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ABSTRACT 


Astronomical data compiled during the last 70 years by the interna- 
tional organizations (ILS-IPMS, BIH) providing the coordinates of the instan- 
taneous pole, clearly shows a persistent drift of the ’’mean pole" {- bary- 
center of the wobble). 

This study was undertaken with a specific objective in mind: to inves- 
tigate the possibility of a true secular motion of the baxycenter; that is, 
an actual movement of the earth pole of figure resulting from differential 
mass displacement due to lithospheric plate rotations. 

The method developed assumes the earth modeled as a mosaic of 
1° X 1® crustal blocks, each one moving independently in accordance with 
their corresponding absolute plate velocities. 

The differential contributions to the earth's second-order tensor of 
inertia were obtained and applied, resulting in nO significant displacement 
of the earth's principal axis. 

In view of the above, the effect that theoretical geophysical models 
for absolute plate velocities may have on an apparent displacement of the 
"mean pole" as a consequence of station drifting was analyzed. 

The investigation also reports new values for the crustal tensor of 
inertia (assuming an ellipsoidal earth) and the orientation of its axis of figure, 
reopening the old speculation of a possible sliding of the whole crust over 
the upper mantle, including the supporting geophysical and astronomic evi- 
dence. 
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"Quand on veut appliquer les theories de la mecanique 
rationnelle a I’etude des phenomenes naturels, on se 
trouve toujonrs en presence de questions d'une compli- 
cation extreme. Si I’on voulait traiter ces questions en 
toute rigueur, on ne pourrait jamais y parvenir, et cela 
pour des raisons diverses que je n’ai pas besoin d'enu- 
merer. Nous sommes done obliges de nous eontenter 
de resoudre, non pas les questions elles-memes que 
nous avons en vue, mais d’autres questions qui s'en 
rapprochent plus ou moins, et qui presentent un degre 
de simplicite assez grand pour que nous puissions en 
aborder la solution plus ou moins rigoureuse. C'est ainsi 
que nous sommes conduits a substituer aux solidea de la 
nature de corps aolides de forme absolument invariable; 
e’est ainsi encore que nous attribuons habituellement aux 
liquides une propie^te de iluidite absolue qui n'existe nulle- 
ment dans la nature, etc. Mais nous ne devons pas perdre 
de vue que, en agissant ainsi, nous nous mettons a cote 
de la realite, et nous devons toujours nous preocuper de 
1' influence que les circonstances dont nous avons fait 
abstraction peuvent avoir sur le resultat auquel nous 
sommes parvenus". 

C. E. Delaunay 

"Sur I’hypothese de la flulditc inter! eure du globe terres- 
tre"; C. R. Acad. Sci. Paris, 67, p. 68 {1868). 
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1 INTRODUCTION 


1. 1 General Background 

The earth, rotation vector is characterized by its magnitude (instantaneous 
spin rate) and its diroction, which can be referred about any prescribed reference 
system. 

Due to the nature of the rotational motions of the earth and other related 
geodynamie phenomena, two basic reference frames are required to fully specify 
the orientation of the instantaneous rotation vector : an exter Jial system with di- 
rections fixed in space (inertial system) and a system rigidly linked to the earth 
(terrestrial system). 

The rigorous definition and realization of these systems is a matter con- 
cerning the International Astronomic Union (lAU), the hiternational Union of 
Geodesy and Geophysics (lUGG) and the International Association of Geodesy 
(lAG). A recent colloquium in Toruh (Poland) was held to initiate and coordinate 
efforts by the international scientific communiiy in the arduous task of defining 
future reliable reference frames (accuracies of 5 x 10”^ radians and about 1 
cm. in orientation and in position, respectively for both the inertial and terres- 
trial systems ) . Diverse viewpoints were presented and debated, but only the 
basic requirements seem to be fairly well established [Koiaczek and Weiffen- 
bach, 1975, p. 13]. 

Since the earth rotation vector continually changes in magnitude and di- 
rection, its study requires the subdivision into three major domains: 

(i) Variations of instantaneous rate of spin- — 

(ii) Inertial orientation 

(iii) Terrestrial orientation 

Table i.l given in [Rochester, 1973] and adapted from [Kaula, 1969] tab- 
ulates the various compononts of the rotational motion of the earth. 
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Table 1. 1 Spectrum of Changes in the Earth's Rotation* 


A. Inertial Orientation of Spin B. Terrestrial Orientation of C. Instantaneous i^in Rate 

Axis Spin Axis (Polar Motion) do about Axis 


1. Steady precession: amplitude 
23®. 5; period s 25,700 years. 

2. Principal nutation; amplitude 

20 (obliquity)} period 18. 6 
years. 

3. Other periodic contributions 
to nutation in obliquily and 
longitude; amplitudes ; 
periods 9. 3 years, annual, 
semiannual, and fortnightly. 

4e. Discrepancy in secular de- 
crease in obliquity; o'\l/ 
century (?). 


1. Secular motion of pole: ir- 
regular, « O". 2 in 70 years. 

2. 'Markowitz' wobble; amplitude 

0*'. 02( ?)} period 24-40 years(7). 

3. Chandler wobble; amplitude 
(variable) s-0''.15; period 
425-440 days; damping time 
10-70 years (?). 

4. Seasonal wobbles: annual, ampli- 
tude “0 ,09; semiannual, ampli- 
tude «o". 01. 

5. Monthly and fortnightly wobbles; 
(theoretical) amplitudes s 0 .'001. 

6 . Nearly diurnal free wobble; 
amplitude < 0" . 02( ?) ; period(s) 
within a few minutes of a sider- 
eal day. 


1. ffecular acceleration; 
u)/t0e;-5 X lO'^Vyr. 

2. Irregular changes: (a) 
over centuries, w/to^ ±5 
xl0”^°/yr; (bj over 

1-10 years, 0 )/co £ ±80 x 
10"’"^yr; (c) over a few 
weeks or months ('abrupt'), 
o)/co^ ±500 X 10 "^^yr. 

3. Short-period variations; 

(a) biermial, amplitude 
»s9msec; (b) annual, 
amplitude s 20-25 msec; 

(c) semiannual, ampli- 
tude 2 * 9 msec; (d) 
monthly and fortnightly, 
amplitudes si msec. 


7. Oppolzer terms; amplitudes s 0^^ . 02; 
periods as for nutations. 


*Erom [Rochester, 1973]. 



The scope of this investigation is restricted, in general, to the terres- 
trial orientation of the instantaneous spin axis (Polar Motion) and in partic- 
ular to its apparent secular motion. 

More than fifteen years have passed since the publication of the mono- 
graph "The Rotation of the Earth: a Geophysical Discussion" by Munk and 
MacDonald [I960]; nevertheless, its contents remain surprisingly up-to- 
date and many of its conclusions substantially valid. The work describes 
and studies diverse geophysical- questions, several still unsolved or without 
completely satisfactory answers. 

The most recent advances in this area were discussed at the second 
Geodesy /Solid-Earth and Ocean Physics (GEOP) research conference on the 
Rotation of the Earth and Polar Motion held at the Ohio State University 
[Mueller, 1975 or Rochester, 1975]. 

Secular motion of the pole is a broad term embracing different pos- 
sible interpretations. It seems necessary at this point to introduce the ex- 
plicit terminology used in this report to avoid any future misinterpretations. 

By secular motion of the pole is commonly understood the displace- 
ment of the astronomically-reduced center of the wobble (mean pole or bary- 
center) referred to the Conventional International Origin (CIO) (see section 
2.4). The interval of time covered by observations is less than 100 years. 

Polar wandering implies a geologic time scale and therefore pos- 
sible large variations in the pole; hence It will not be the primary concern 
of this work. 

The astronomically-derived secular fluctuations of the pole, as is 
well recognized, may in fact be 

(^) true secular motion 

(b) apparent secular motion caused by systenmmtic-drifting of the 
observing stations or variations in their deflections of vertical 

(c) true and apparent effects combined 
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There are speculations about the possibility of a true secular dis- 
placement of the pole, at this time stiU a controversial question, with sev- 
eral well-known authors for and against it. Observation data clearly shows 
a persistent drift of the barycenter, but the reasons for this drift are not 
yet clear. Recent geophysical theories have been postulated to explain this 
phenomenon. In the following two sections the available astronomical evi- 
dence and different geophysical interpretations wOi be reviewed. 

1.2 Astronomical Evidence 


Although the possibility of displacement of the earth pole of figure 
due to geological causes (earth mass displacements) was recognized a long 
time ago by geophysicists [Darwin, 1877], geodesists [Helmert, 1880, p. 420] 
and astronomers [Tisserand, 1891, p. 485], the observational verification 
needed to wait xmtil accurate instruments and methods were available. 

Several erroneous rates for the secular polar motion were reported 
prematurely, probably caused bj over-anxiety to Second the theoreticians. 

The indisputable fact that the pole has a progressive motion was established 
after the founding of the International Latitude Service (ILS) in 1899. The 
ILS uses five stations located at nearly the same latitude (38® 8^N) allow- 
ing them to observe the same stars, thus eliminating the effects of certain, 
systematic errors in the star catalogs [MueJier, 1969, p. 81]. 

Lambert [1922] was the first to report a definite polar drift derived 
from ILS data. This initial estimate marks the beginning of an international 
campaign aimed at detecting the reality of a secular displacement of the pole. 

Table 1.2, adapted from [Poma and Proverbio, 1976], presents various 
rates and directions of the secular drift of the mean pole based on ILS ob- 
servations- 


Table 1. 2 


Values of the Secular Motion and Direction of the Mean Pole 
Determined by Different Authors Based on ILS Latitude Observations* 


Author 

Year 

Angular 

Velocity 

Direction 
(Longitude W) 

Period of 
Observations 

Lambert, W, D. 

1922 

0".0066/yr 

83° 

1900-1918 

Kimura, H. 

1924 

58 

67 

1900-1924 

Wanach, B. 

1927 

47 

42 

1900-1925 

Wanach, B. and Mai'inkopf, N. 

1932 

63 

84 

1900-19U 

Wanach, B. and Mahhkopf , N, 

1932 

51 

62 

1900-1923 

Hattori, T» 

1947 

45 

73 

1900-1939 

Hattori, T. 

1959 

359 

65.2 

1900-1947 

Markowitz, W. 

1960 

32 

60 

1900-1959 

Yumi, S. andWako, Y. 

1966 

269 

72.9 

1933-1966 

Stoyko, A. 

1967 

32 

69.8 

1890-1966 

Proverbio, E. et al. 

1971 

294 

65.6 

1900-1962 

Proverbio, E. etal. 

1972 

277 

60.7 

1900-1962 

Proverbio, E. andQuesada, V. 

1973 

307 

69.6 

1900-1969 


*From [Poma and Pr overbio, 1976] 



After reviewing all the available polar motion data, Sekiguehi [1954] 
showed for the first time that besides a progressive drift of the pole, some 
random discontinuities were also present, hi 1960 when the secular motion 
of the barycenter seemed to be well confirmed, an empirically derived li- 
brational component (about a 24-year period) was reported by Markovita 
[I960]. 

The ILS was reorganized in 1962 into the International Polar Motion 
Service (IPMS). Since then it provides the fundamental pole path derived 
from the five initial stations, and also publishes a pole path based on varia- 
tions of latitude received from about 40 other cooperating stations (see 
Table 4.10)that observe different stars [Yumi, 1975]. 

Since 1955 another service, the Bureau International de I'Heure (BIH) 
has also reduced latitude observations. Since 1967 measurements of time 
and latitude from about 50 observatories have been used simultaneously for 
computing the coordinates of the pole, while previously only latitude deter- 
minations were utilized. Unlike the IPMS, the BIH deals With problems con- 
cerning the preservation of its newly-defined system (1968 BIH system), 
therefore providing the possibility of adding or cancelling cooperating obser- 
vatories [Guinot et al. , 1970] . 

A recent stu^ by [Poma and Proverbio, 1976] shows tmquestionably 
that the secular drift of the barycenter, as deduced from BIH polar coordi- 
nates, corresponds completely to that obtained with IPMS data, corroborat- 
ing the reality of a secular motion of the mean pole. However, as mentioned 
above, because of crustal motions, etc. , the observed secular drift of the 
pole may not represent a true polar drift. 

1-3 Geophysical Theories 

The intuitive assertion that changes in the earth’s mass distribution 
may produce displacements in the direction of the earth's principal axes was 
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originally formulated analytically by Darwin [1877], The failure to detect 
any sensible polar shift through astronomical observations diminished the 

'I 

interest in the subject vmtil 1922, when the first estimates of secular drift 
of the mean pole were reported. 

Later, around 1950, the question was very much reopened after devel- 
opments in paleo magnetism suggested that large-scale displacements of 
the rotation pole had occurred. Extensive material has been published show- 
ing different pole positions according to the location of the magnetic rock 
specimens. A recent investigation suggests that the polar wander path as 
derived from paleomagnetic data reilocts mainly continental motions and 
not true polar shift, which remains within 4“ in 55 m.y. [Jurdy and Van 
der Voo, 1974]. 

The first modern considerations to explain a possible drift of the 
pole are due to Gold [1955], who following Darwin's reasoning, postulates 
the great importance of mass redistributions in the changes of the earth 
pole of figure. 

Burgers [1955] in his analytical treatment, confirmed Gold's conclu- 
sions, primarily that asymmetric contributions to the moments and products 
of inertia of the earth are the main motive for polar wandering. 

Several papers including [Pekeris, 1935] and [Runcorn, 1962] have 
suggested mantle convection as a possible cause for polar wandering. 

Munk and MacDonald [I960] extensively reviewed the pros and cons 
of polar wandering, using different theoretical models. After claiming that 
fcom examination of various possible excitations of polar wandering it ap- 
pears that the distribution of continents and oceans exceeds other possible 
causes by a factor of 1,000, their final conclusions recognize the complex- 
ity of the problem (polar wandering), stating that it remains unsolved. 

New support to the theory that Inertia changes in the earth may pro- 
duce large-scale polar wandering was Introduced in the often-quoted paper 
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by Goldreich and Toomre p.969]. 

TakeucM and iSugi [1972], using a crude calculation based on Darwin’s 
equations (see section 4.5), concluded that mantle conveetioh may cause 
secular changes in the earth products of inertia and consequently polar dis- 
placement. The direction of wandering calculated is about 90° W but its 
rate is about three times as large as the astronomically observed value. 

While this work was in progress [Liu et al. , 1974], using Darwin's 
approach, obtained changes in the products of inertia of the earth due to 
new creation of lithospheric material and found a secular shift of the pole 
toward the direction 67° W although only 10% of the astronomically observed 
rate. 

Finally [Pan, 1975] speculates that polar wandering is dynamically 
linked with the secular tectonic movements in the earth’s upper layers. 

Thus, as can be seen, the majority of authors tends to agree with 
one hypothesis; changes in the inertia tensor of the earth may produce secu- 
lar motion of the pole on a short term basis and polar wandering, if large 
redistribution of masses occurs. 

1.4 Scope of the Ihvestigation 

The main objective of the present study is to calculate the effect of 
crustal mass displacements on the secular motion of the earth’s pole within 
the past 70 years. 

The differential contributions to the earth tensor of inertia from 
differential tectonic plate rotations are evaluated to verify the geophysical 
theories . 

Borrowing the following quotation from Jeffreys [1974], 

"A theoiy consists of three parts; 

(1) a hypothesis p is set up for consideration; 


(2) the consequences of p are worked out, as for example, qi,qs. • . 

(3) finally, qi,qs. . . are compared directly with observations. 

Part (2) is generally mathematical; part (3) may involve considerable 
statistics. 

If all qi ,q 2 . . . are found to agree with observations, we may say 
that p is strongly confirmed. If some agree but existing observa- 
tions do not test the others, there is some confirmation, but less 
strong. If even one is definitely contradicted by observation, p 
must be either rejected or modified”. 

This investigation is concerned primarily with point (2). The hypo- 
thesis considered is the inlluence that tectonic plate motions may have on 
the true or apparent secular motion of the polo. 

At this point one should consider that the astronomical observations 
cover an interval of time of less than a century, while the absolute plate 
velocity models used in this investigation have implicit longer-term rates. 
Thus the observations essentially reflect possible current changes of tho mean 
pole, while the geophysical plate models involve velocity rates given for 
one year but obtained from data covering a geological time span. 

The following chapter studies some basic geodynamic problems after 
reviewing the rigid earth case. The necessary theory for evaluating the 
differential contributions to the earth tensor of inoi'tia due to differential 
mass displacements is introduced, A complete formulation for transform- 
ing tensors of inertia is included in Appendix B, where a novel approach for 
matrix diagonalization is presented. Contrary to the usual procedures, the 
orthogonal matrix of the transformation is obtained first, from which the 
eigenvalues are easily computed. Other derivations relevant to this as well 
as the other chapters are given in the different appendices in a convenient 
unified matrix notation. 

The third chapter describes the various models adopted in the inves- 
tigation. These include a crustal model based on the theory of isostatie 
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compensation, plus tlie assumed break-up of tectonic plates and their abso- 
lute velocities with respect to the underlying mantle. 

Chapter 4 is the presentation of the numerical experiments and the 
results obtained. The effect on the pole displacement of changes in inertia 
or possible station drifting is analyzed, . 

The last chapter contains the summary of the investigation and the 
conclusions reached. 


2 GEODYNAMICS 

The study of the rotation of the earth about its center of mass has 
constituted a fundamental problem in astronomy for centuries. Recent im- 
provements in the theories of earth mechanics and the development of instru- 
ments and new observing techniques have forced the inclusion of the topic 
in the domains of geophysics. 

Unfortunately the subject matter is greatly complicated by the het- 
erogeneity and deformabiliiy of our planet, plus the action of changing forces 
mainly due to the attraction of other bodies in the solar system. 

The effect on the rotational dynamics of the earth of such varied 
phenomena as winds and moving air masses; atmospheric, oceanic and bod- 
ily tides; sea level changes and tectonic motions; rigidity, elasticity or 
plasticity of the earth mantle, the nature of its core and their interactions, 
are not yet fully understood. The tact is, earth dynamics remains a clialleng- 
ing area for investigation proved by the intensive world-wide theoretical and 
applied research efforts devoted to the subject. 

This chapter wiU present alternatives for formulating the general 
Lagrange- Liouville equations after consideration of specific geodynamic 
hypotheses. Some known kinematical and dynamical postulates will be re- 
viewed for use in subsequent chapters. 

For its simplicily and great didactic importance, the ideal case of 
rigid earth will be treated first. Besides the necessary formulation, this 
should provide a conceptual guidance for understanding ulterior developments. 
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2.1 Preliminaxy Coccepts and Notation 


To facilitate writing mathematical expressions, a simplified matrix 
notation will be used throughout this work. In this way the explicit represen- 
tation of aU elements in the equations is avoided. 

The benefits of matrix calculus over vector or tensor methods are 
obvious in a computationally-oriented scientific community. The resulting 
formulas can be programmed immediately, taking advantage of matrix alge- 
bra and the subroutines available in most computer centers . 

Unless otherwise stated, the following notation will be used in the 
present and remaining chapters: 

2. 1. 1 Coordinate J^stems 

Right-handed rectangular coordinate systems will be represented in 
general by the convention 

(Xi,X2,X3) = (x) (2.1-1) 

with the implied origin denoted by Ox. In dynamics , when the system has 
the center of mass (CM) of a body as its origin, it will be called a "central 
system" and the following notation applies 

(XoXjXo2,Xo 3) = (Xo) 

2. 1. 2 Vectors 

In the treatment of vectors, several approaches are possible. In 
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this study vectors referred to an (x) coordinate system will be presented 
by column matrices denoted as 


[a}. = 


ai 

a.2 





( 2 . 1 - 2 ) 


A vector of this lype is called a column vector and defines a point in the 
Euclidean space which can be specified by three numbers a, {i = 1,2,3} 
expressing the coordinates along the Cartesian frame (x). 

To conform with matrix multiplication rules, sometimes the compon- 
ents of the same vector {a}^ will be arranged in horizontal array 


{a]J= [ai as as}^ (2.1-3) 

and termed row vector. 

One-dimensional vectors are called scalars. 

In cases where the coordinate system used is clear by the context, 
the siibindex x in the vectors will be omitted. 


2.1.3 Matrices 


In general a 3 X 3 real matrix will be denoted symbolically by f M] . 
Nevertheless throughout the text some well-lmown types of matrices are writ- 
ten without the brackets. For example, this is the case of the rotation matrix 
R between two Cartesian systems of coordinates. 

The following special types of matrix notations are introduced and 
subsequently used: 

a) Skew-symmetric (anti-symmetric) 

To every vector [a]^ it is possible to associate a skew-symmetric 
matrix denoted by 
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where clearly; 


[all = -[ak (2.1-5) 

The matrix in (2.1-4) is also referred to in mathematical literature as "axial 
vector" or "pseudovector". Anti-symmetric transformations are reviewed 
in Appendix A. 

b) Identiiy matrix 

The 3X3 unit matrix will always be denoted by [1]. 

c) ^mmetrie 

Such matrices usually will be denoted by the upper triangular elements 
only, and an "s" in the lower left of the matrix. 

d) Diagonal 

The symbolic notation [iVLJ will apply to diagonal matrices excluding 
the unit matrix. 


2. 1-4 The Tensor-Matrix Concept 


hi general, an operator , such that 
lCO 


(X) 


-> (x') 


( 2 . 1 - 6 ) 


is called linear if V{x} , {y} £ (x) and any real nuniber k, the following 
two properties are satisfied , 

= {x} (2.1-7) 

and 


L^’^^Cx+y] = L ^ Cx} + {y} 


= T 


(r) 


( 2 . 1 - 8 ) 



A linear operator L ^ { L r scalar) is called a tensor of order 

(rank) r if 

V r >1 and V {xj 

the application L ^ {x} is a tensor of rank r - 1. 

A consequence of the above definition is that 

(i) Any vector is a tensor of order 1, since when applied to other 
vectors it gives a scalar, i. e. , 

[y} = scalar (order = 0) [x] is a tensor of order 1 . 

(ii) Matrices are tensors of order 2, that is, 

[ M ] [x] = vector (order = 1) =- • [ M ] is a tensor of order 2. 

Accordingly, 3X3 matrices obtained through the components of any 
vector [x] such as {x} [x}’^ and [jc] [x]^ are tensors of second order. 

The term "dyadic, " infrequently used today, occasionally is mentioned synon- 
imously with the {x}[x}^ tensor. 

In general the number of components in a tensor of order r and dimen- 
sion n is given by n*'. This study will be limited to tensors of r £ 2 and n = 3. 

At this point one should notice the diffex'ence between rank of a tensor 
and ranic of a matrix (Rank [ M ] = the number of linearly independent columns 
of [M]. Obviously 1 ^ Rank [M] < 3). 

As a genei’al rule, the algebra that has boon developed for matrices 
may be used for tensors as well. In particular, tensors of second order, 
like matrices , maybe diagonal, symmetric, anti-symmetric, orthogonal, etc. 

Finally, the process by which from a tensor of order r ^ 2, other 
tensors of order r - 2 are derived, is Icnown as "contraction" of a tensor. 

The contraction of a tensor of order 2 is a scalar equivalent to the trace of 
its matrix. 


2. 2 Inertia Tensors 


There is little agreement in physical and mathematical literature 
about what is referred to as the tensor of inertia. Therefore, in order 
to avoid confusion, the foUowir^ definitions and notations will be used: 

2 . 2 . 1 Associate Tensor of Inertia (Associate Inertia Matrix) 


The associate tensor of inertia of a body with respect to a system 
(x) with origin Ox is defined by 


r 


-X1X3 

-X1X3 


A 

-P -E 

[I]= / IxKxrdm = | 


Xi+: 3 ^ 

-XgXa 

dm = 


B -D 


f 

a 

a 

xf+4_ 


3 

— 

C 


( 2 . 2 - 1 ) 


where A, B and C are referred to as the “moments of inertia" of the body 
with respect to the xi , X 3 and Xg axes respectively, and D, E and F are 
called "products of inertia" with regard to the axes Xg & Xg , Xi & Xg and 
Xi & X 3 respectively. The set of moments and products of inertia are known 
as the six constants of the bo(^ with respect to the (x) system. 

Considering ( 2 . 1-5) one may write: 

[ [x][xr Y = Ul[x]^ = [x]^[xl = “[x][x] = -[xl® (2.2-2) 

Thus the associate inertia tensor is clearly symmetric. 

By an orthogonal transformation (which preserves the rectangular 
Cartesian character of the coordinates), the matrix [ I ] may be reduced to 
the diagonal form 
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Ap 0 0 

'Bp 0 

5 Cp 


(2.2-3) 


Then the ’?alues Ap, Bp, Cp are called the '’principal moments of inertia" 
(with respect to the reference point OxsOXp) and the (Xp) axes of the trans- 
formed coordinate system are called the "principal axes of inertia". 

The tensor of inertia referred to a (Xq) system with origin at the 
CM will be called "central tensor of inertia". If the central tensor is prin- 
cipal (i. e., diagonal) its associate system is called "central principal axes" 
or sometimes "axes of figure" of the body. They will be represented by 
the notation 

(Xopl » XopS» ^0p3 ) = (Xop) 

Observe that according to the above definitions the principal axes of 
inertia of a body are not necessarily central. 


2. 2.2 The Tensor of Inertia (Inertia Matrix) 


The tensor of inertia of a body with respect to a system (x) with 
origin Ox is defined as 


m 



(2.2-4) 

r E 

Egg D 

E33 


where In, Hag, Has are referred Lo as the moments of inertia with regard 
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to the planes Xi = 0, Xs = 0 and xa = 0 respectively, and D, E and P the 
products of inertia with respect to planes Xs = 0& Xa=0, xi = 0& X 3==0 
and Xi = 0 & Xa = 0 . 


The symmetry of [ I ] follows immediately from the matrix equation 

[tx}Cx]^]^ = {x}{xr 

Finally the moment of inertia with respect to the origin of the coordi- 
nate system is given by the contraction of the tensor [ I ] or 


tr[lE] 


= fixV{x]dxa = Hxi+ xl+ 3!^)dm 


(2.2-5) 


In Appendix B the reader may consult the different transformations 
between the tensors [ I ] and [ 1 1 as well as the translation and rotational 
effects on the tensors of inertia. 


2.3 Rigid Earth Case 

Assumption 1: The earth is rigid. 

Assumption 2: The central principal axes of the rigid earth (Xop) 
are taken as system of reference. This system is body-fbced but moves with 
the earth in space. Under the above two hypotheses, the Lagrange-Liouville 
equations given in Appendix D by the general expression (D.2-12) will reduce 
to the well-known Euler's dynamical equations. 

The following simplifications are a consequence of the two assump- 
tions established; 

1 [io] = [0] 

[ AI] = [0] 

{6£o} = {O} 

Ch3 = {o}==> [h] = [ 0 ] 

{Ls} = [0] 

A 2) Reference system; (xop) ==> Do = Eq = Po = 0 [Iq 1 [“Iq-J 


I 


] 

j 


i 




I 
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(2.3-1) 


Hence, equation (D.2-12) reduces to 

{L} = t'IoJ{cl33 -i-[a3]Mo.l{w] 

After the matrix multiplications are performed, the above equation may be 
explicitly written as 


= Acci 

- (B-C) 0 ) 30)3 

(2.3-2a) 

= BW 3 

- (C-A) 0)3 0 )t 

(2.3-2b) 

= C CO 3 

- (A-B) 0 )i 0)3 

(2.3-2C) 


where all the vector components and the moments of inei'tia clearly refer to 
the selected reference system (xop). 

2.3.1 Point of View from the Earth Surface 

A variant of equations (2.3-2), interesting from the point of view of 
an observer on the earth, is obtained under the following assumption: 

Assumption 3: The body torques acting on the eiirth are scero. 

The solution of the differential equations (2.3-2) when {L) = {O] 
is well documented in the classic works on dynamics (French literature 
Gonuuonly uses the notation p - q and r W 3 ). 

An analytic solution of ( 2 . 3-2) under the constraint imposed by the 
third assumption is possible. The obtained values of {oJ } ai‘e expressions 
containing Jacobi’s elliptic functions. For a detailed discussion and exten- 
sive bibliography on the subject, consult [Leimanis, 1965]. 

Nevertheless, similar final conclusions may be derived from the 
simplified case when two principal moments of inertm are equal. This is 
a rather good approximation since the earth is nearly axially symmetric. 
Thus, 

Assumption 4: A = B. 
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The above assumption implies that all the lines in the principal equa- 
torial plane (plane through the CM normal to the axis of figure, also referred 
to as "dynamical" equatorial plane) are principal. Therefore the first two 
axes of the reference system are not uniquely defined. In the analysis that 
follows the original Cartesian frame (xop) will be kept as reference system. 
That is, the axes of reference are the principal axes of the rigid earth cen- 
tral momenta! ellipsoid satisfying the condition A < B < C. 

Thus finally, equations (2.3-2), after assumptions 3 and 4 are imple- 
mented, become 

+ ncog = 0 (2.3-3a) 

0)2 - nWi = 0 (2.3-3b) 

0)3 ~ 0 (2.3— 3c) 

where the following substitution has been introduced 
C-A 

n = (2.3-4) 


Equation (2.3-3c) immediately gives 

0)3 (t) = 0)3 (0) = 0)3 H constant (2. 3-5) 

Changing to complex variable notation, and defining 

0 )o = Oil + io)a (2.3-6) 

the remaining two differential equations (2.3-3a & b) catibs combined in a 
single differential expression in O),, 

0 )e - inoJo = 0 (2,3-7) 


Equation (2. 3-7) has the exponential solution 
... i(nt+8) 


(2. 3-8) 
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where a and jS are two constants of integration. 

Equation (2.3-8) yields the known periodic solutions: 

ct>i(t) = acos(nt + j8) (2. 3-9a) 

C0a(t) = o'sin(nt+ J3) (2.3-9b) 

The earth rate of spin is the magnitude of the angular velocity vector [o)}. 

|io| = V == V 4- COa constant (2.3-10) 

From (2. 3-6) and (2. 3-9) one may conclude that the vector oOo describes 
a circle of radius CL in the principal equatorial plane while 0 O 3 remains con- 
stant. Hence {to] precesses uniformly in the body frame about the symmetiy 
axis with angular velocity n (see Pigs. 2. 1 and 2.2). Once the differential 
equations are solved, other important quantities related to the rotation vec- 
tor, such as its direction cosines, can be obtained immediately (see Fig. 2.1). 

tOi 

mi = — r (i = 1,2,3) (2.3-11) 

I wl 

The coordinates of the instantaneous rotation axis Pj in a plane paraJUol 
to the dynamical equator tlirough the north "inertia pole" of the earth (Inertia 
pole = point of intersection of the earth with its axis of figure), and referred 
to the axes of reference Xopi and Xopa are 

COi 

Xopi = Up (2,3-12a) 

CO3 

CJs 

^op3 ~ (2.3— 12b) 

CO 3 

where 

Tp semiminor earth axis 

The rate at which the angular velocity vector (instantaneous spin axis) turns 
about Live Xop 3 axis was given by 
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Xbpi(A) 


Xbp 3=0 3irp=prijicipal equatorial plane 


Fig. 2.1. Bistantaueous Angular Velocity Vector {to} in. 
the Central Principal Frame (xop) 






(2. 3-13) 


The moments of inertia of the earth are approximately related by 

C“ A - rtrinn 
— - — fa 0. 0033 
A 

The component CO 3 of the angular velocity is for all practical purposes equal 
to ( coj (one rotation per day); thus finally, 
n fa 0.0033 rotations per day 

Therefore, the period of {to} around the rigid earth axis of figure is 

T = - » 305 sidereal days (2. 3-14) 

n 

This is generally termed the Eulerian period. In conclusion the following 
dynamical postulate may be stated when assumptions 1 to 4 arc fulfilled: 

Postulate 1 : An observer at rest on the earth will see that the instan- 
taneous rotation axis intersects its surface tracing out a circle (polhodc) 
around the reference pole (mean pole ^ pole of figure) and moving in a retro- 
grade direction with a period of about 305 sidereal days. 

The polhode may also be referred to as wobble or free wobble. 

As a corollary, one has: The mean position (center*) of the wobble 
described by the spin axis coincides^ rigorously with the pole of figure. 

Sometimes the name "free nutation" is applied to this motion on 
the earth of the instantaneous angular velocity vector. This denomination 
is nevertheless unfortunate, as rightly pointed out by Rochester et al. , 

[1974] and should be applied only in connection with a space-fixed system, 
as is implied by the term "nutation". 

The fact that the earth is not rigid introduces further difficulties 
discussed later in this chapter. 


2.3.2 Point of View from Space 


Up to now the described motion of the instantaneous spin axis of the 
rigid earth conforms with an observer at rest in the rotating body-fixed 
reference frame. 

The situation changes if one considers an external observer related 
to an inertial (space-fixed) system - 

As is known (see section D. 1) an immediate consequence of the assump- 
tion {L j = { 0 } is that the total angular momentum vector [h] is constant 
in the inertial frame. Hence the direction of {H} will remain truly fixed 
in space and for this reason is also called the "invariable axis". The plane 
through the earth CM normal to [H] is termed the "invariable plane". 

The relation between [H} and the earth angular rotation vector [co} 

(the components of this vector are known for any particular instant at which 
the differential equations of motion are solved) is given by equation (D. 1-22) , 
namely: 

{H3 = [Iol{to}+ {h} (2.3-15) 

Thus the direction cosines of the vector [H} with respect to the earth-fixed 
system are immediately available: 

hi = — (i = 1,2,3) (2.3-16) 

iHl 

After consideration of assumptions 1 and 2 equation (2. 3-15) reduces 


to 


Ch3 = rioJCco] 

(2.3-lY) 

or 


Hx = Aoh 

(2. 3- 18a) 


Hs = BtOs 

(2. 3-18b) 
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Hs “ C CiJs 


(2.3-18C) 


and 

= A®co!+ B®oo|+ Cfwl 


(2.3-18d) 


Denoting by ( 0, 0, 1) the components of a unit vector along the earth axis 
of figure, the condition of coplanariiy between this principal line and the 
vectors [co] and {h} is 


0 

COi 

Acoi 


0 

oiz 


1 

t*J3 

CdOa 


= 0 


(2.3-19) 


or 


BtOi COa - AcOi COg = 0 


(2. 3-20) 


The above equation implies that when the fourth assumption (A - B) is satis- 
fied, the three axes (rotation, figure and angular momentum) are in the same 
plane at any instant. Since [h] is constant, the plane rotates around the 
invariable axis. 

The axial vector representing the angular momentum forms an angle 
with the earth axis of figure defined by (see Figs. 2.3 and 2.4) 

Vh?+ Hi 

tan0 •= (2.3-21) 

H3 

where the components of [H] in general may be computed from equation 
(2.3-15). 

If rigidity and A = B are assumed, one has 
VA^(tOi+ o)|) Aa 

tan0 = - — " constant (2.3-22) 

CtOa C6U3 


Hence the following postulate can be stated; 

Postulate 2; The axis of figure describes a cone in space whose axis 
is the invarUtblo line ([ll] ) and whoso somiangle is 0. 
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The angle £ between the instantaneous spin axis and the earth axis 
if figure is given by the relation 


tan ^ = 




(4 


60s 


(2.3-23) 


and under the four mentioned assumptions, it follows 


tan C - — = constant 

60 s 


(2.3-24) 


Thus the instantaneous rotation axis forms a cone (body cone) of 
constant angle ^ with the axis of figure. This was also implied previously 
by equations (2. 3-5) and (2. 3-9) and was summarized, in Postulate 1. 

As a consequence, one has the following corollary: The invariable 
line describes a cone in the body whose axis is the axis of figure and whose 
semiangle is 0 . 

Prom (2. 3-22) and (2. 3-24) it is obvious that for a rigid earth and 

A = B, ^ 

tan C = T 0 
^ A 


The above implies that the displacement of the pole of figure from the pole 
of rotation is greater than its displacement from the pole of angular momentum 

Q 

by a factor r « 1. 00327. 

i A 

Finally, the angle between the angular momentum vector and the 
instantaneous rotation axis may be computed through the equation 


cos (A HOP| ) 


{Hrcco} 

([Hr CH})^®({cor [to})^® 


(2. 3-25) 


In the particular case of rigidity and A = B 

Aa^+ Cool 

cos (A HOP] ) = -- - — ' = constant (2. 3-26) 

Va®q?+C®co 1 ^ O'® + 6o| 
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Postulate 3 : The instantaneous axis describes a cone in space (space 
cone) whose axis is the fixed angular raomentum axis and whose semiangle 

is e - e. 

This motion of the instantaneous earth rotation axis around the invar- 
iable axis, even in the absence of torques, has been called "sway” by Munk 
and MacDonald [1960, p. 45], but as opportunely mentioned [Rochester et 
al., 1975] the name "free nutation" is mor.' appropriate. 

Clearly if A = B = C, 

cos [L HOP, ) = !==>/- HOP, = 0 
This implies the vectors [h] and [io} have the same direction. 

Therefore, the instantaneous rotation axis and the total angular mo- 
mentum axes are coincident only when A = B = C. (The momental ellipsoid 
is a sphere and every diameter is a principal axis). 

Moreover, one has 

A=B=C=^ n =0 

Thus the spin axis is fixed in space and in the body. 

Equations similar to (2. 3-25) may be applied to obtain the cosine of 
the angles 6 and £ 

Cto,'^ 

cos 6 = — — (2.3-27a) 

CO 3 

cos ^ = - - - - (2.3-27b) 

■\ja^+ co| 

Therefore, recalling ( 2 . 3 - 22 ) and (2.3-24) 

OS 

sin 0 = ■ — - ■ (2.3-28a) 

\IA® 0 '®+ C^co| 

sittC = - - — (2.3-28b) 

V t- tOfl 
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Thus, finally one can write (2.3-26) as 

coa C CO3 a 

cos (A HOP,) = - 3 ^ - - + ^ 

^ja^+0l3% sja^+icii 

- cos Ceos 6+ sin C sin 0 = cos{C - S) 
and 

Z.HOP, = C - 0 


Aa 

■v/A®a®+C®col 


(2.3-29) 

(2.3-30) 


Similarly, 
sin (A HOP,) = 


When the angles are small, as in this case, 

sin(AHOIj) « AHOI] , cos0s=rfl , sinC«C 
From astronomic observations. 


sin{ C “ 6) = sin C cos 0 - cos C sin 0 

C cOa C— A C— A 

= cos0sinC (2.3-31) 

Va^q:®+ C®co| -^a^+col C C 


H = ~ 0,003272 = dynamical ellipti city 

and 

C < 0".3 fa 9.5 m 

Thus 

A HOP, ja 0 ". 00098 fa 3.108 cm 


(2. 3-32) 

(2.3-33) 
(2. 3-34) 


Let an inertial system (X) be introduced with origin OX = CM and the third 
axis in the direction of the fixed angular momentum vector {h} , This sys- 
tem can be related to the earth-fixed central principal axes through the 
Euler angles <p, 0, t/> (see Fig, 2.3) through the transformation 

tx}=]a[xop) (2.3-35) 


where the rotation matrix IR is given by 
H= R3(-tp)Ri(~e)R3(~il)) 
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(2.3-36) 


or 

Ixop} = (2.a-37) 

where 

= R3{0)JRi(0)R3{(p) 

In this case, the Euler kinematic equations are given by (see Appendix C) 


tol 


sin0 sin ij) 

cosij} 

0 * 

Wg j 

> = 

sin0 coscf) 

-simp 

0 

1 

tOa 

1 ‘ 1 


cosd 

0 

1 


(2.3-38) 


* ■ » 

Clearly, the values of the angular velocities cp, G and i/) will give the rate 
of precession, nutation and spin of the earth-fixed system (xop ) with respect 
to the inertial frame (X). 

Applying the above to the ideal case of rigidity and A = B, it follows 
immediately from (2. 3-22) tliat 0=0. 

Substituting the values of (tOi,cOa,cOa) from (2.3-5) and (2.3-9) in 


(2.3-38), 

<Osin0sini/) = asin(nt>^) 

k 

(psinGcosi^ = 0£cos(nt + ^) 

« « 

l/) + cp cos 0 = C0.3 

From the first two equations it follows immediately that 

% # 

(p^in^ 0 = ==> cp sin 6 = ± a 


(2- 3-39a) 
(2.3-39b) 


(2.3-39C) 


(2.3-40) 


Substituting the above in (2. 3-39a) one has 

cjj = nt -i-j8 ^ = n (the rate of spin) 


(2.3-41) 


where n is given by (2.3-4). 

Replacing (2.3-41) in (2.3-39c) 


cp cos 0 = CO3 - n 



C-A 

A 


) 



(2.3-42) 
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Fig. 2.3. Euler's Angles between the Inertial (X) 
and Central Principal (xbp) Frames 



Fig. 2.4. Earth's Spatial Motion 


and therefore, using (2. 3-40) and (2. 3-42) 

cp^ = (2.3-43) 

and 

fV 

tan 0 = — — = constant 

C 0 O 3 

as was found before in (2. 3-22). 

The result of equation (2,3-41) may be interpreted kinematically as 
the rolling of the body cone on the fixed-space cone without slipping (see Fig. 
2 . 3) at a rate of spin n precisely the same as Ibund before for co when it 
was referred to an earth-fixed frame [Poinsot, 1851]. 

In order to visualize the complex behavior of the motion of the earth 
even under the simple assumption of rigidity, the description of Fig. 2.4 
may be illustrative. 

The position of the body cone has been drawn at two different epochs 
Ti (initial epoch) 

Ts = Ti + J (T~ 305 days) 

The invariable axis remains fixed in space at all times when [L) = {O} . 

Attached to the body cone is represented an ideal observatory, the 
position of which with respect to the rigid earth will not change from epochs 
Ti to T 3 . Nevertheless, it will wobble in space as the earth does. (Notice 
that for the purpose of clarity the figure is not drawn to scale ). 

The graph shows clearly the variation of observed instantaneous co- 
latitude at epochs Ti and Tg , while the reduced colatitude (referred to the 
axis of figure or mean pole) remains constant. Thus knowing the position 
of the instantaneous rotation axis with respect to the axis of figurii (polar 
motion), one may obtain the observatory reduced colatitude- (In practice 
a conventional "geographic” referoneo system is used, although the same 
conclusions apply. In this case reduced colatitude will apply to ilie third 
axis selected frame). 
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Note that when two observatories are 180*’ apart, the observed cola- 
titude variations should be opposite in phase, that is, at the same epoch 
Ti , maximum colatitude in one observatoiy will correspond to minimtim 
Golatitude in the other. As it turned out, this was the way polar motion 
was corroborated in 1891 after observations in Hawaii and Berlin were pro- 
cessed. For a complete account of the intricate history of the discovery of 
latitude variations, the review in [Lambert et al. , 19311 is recommended. 

The methods currently used for obtaining latitude are properly explained in 
[Mueller, 1969] and [Melchior, 1957] where many other references on the 
subject are available. For specific data reduction procedures at the appro- 
priate international organizations, consult the annual reports of the BIH 
and of the IPMS. 

Up to now the condition {l} = {O] was enforced. It is quite well 
known that this is not the ease and that the action of external body torques 
(primarily due to the gravitational effects of the sun, moon and planets) cause 
motions of the invariable axis in the inertial system. Consequently the vec- 
tor H will process and nutate in space. Other periodic contributions to nu- 
tation (i. e. , annual, semiannual) are possible because the external body 
forces acting on the earth depend on the ephemeris of the celestial bodies 
where clearly several periodicities are present. 

No further considerations will be devoted to the ease when the body 
torques are not zero. The interested reader may consult classical works 
like [Usserand, 1870] or the modern treatment by Woolard [1953]. 

2.4 Non-Higid Earth 

Although in the previous section the hypothesis that the earth is a 
rigid body was considered, some given functional relations using [h] or 
{(o} are general and may be applicable to any circumstance in which these 
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vectors are known. 

Under the assumption of absolute rigidily, the relative positions of 
all mass particles are constant, so the external form as well as the value 
of the central (principal) tensor of inertia is fixed and independent of time. 

It has become increasingly clear that the assumption of rigidity for 
the earth is incorrect and a better modeling should be investigated. Because 
the mass distribution of our planet is subject to variations with time (tidal 
deformation, crustal motions, etc.) producing changes in its central inertia 
tensor, the rotational dynamics of the earth is better studied by the general 
Lagrange-Liouville equations as given by (D.2-9), (D.2-12) or (D. 2-15), 
depending on the reference system adopted. In the case of a deformable body, 
the total action of the torques is expressed by the sum [L] i- [Lg } where 
[L] are the body torques and [Lg} the surface torques (see Appendix D). 

Theoretically an infinite number of reference systems are possible, 
although it is convenient to select the one best suited to its specific purpose 
and which is practically feasible. 

Previously a theoretical reference system was selected, the one de- 
fined by the central principal ivxes of the rigid earth. Considering that this 
frame is not available through observation, a crust-fixed or ''geographic” 
system is generally adopted. 

The present choice of this type of a terrestrial system conventionaUv 
used by geophysisists , astronomers and geodesists, is defined as follows 
(see Fig. 2.5); 

Origin: Close as possible to the geocenter (center of mass of the earth in- 
cluding the atmosphere). 

Xg axis; Directed toward the CIO (Conventional international Origin) as 

defined by the IPMS (International Polar Motion Service) and the 
BUT (Bureau International do I'lleuro). Some doubts about the defi- 
nition of the CIO had been raised and other {lossible redefiniUons 

\V.\ 




mentioned [Melcbior, 1975]. 

Xi axis: Passes through the point of zero longitude as defined by the 1968 
BIH system [Guinot et al. , 1971], 

xg axis: Forms a right-handed coordinate system with the xi and Xa axes. 

Clearly, for short periods of time ('-> 1 year) the above conventional 
system may be considered earth-fixed. H: is known that rigorously this Is 
not the ease, due primarily to crustal tectonic motions which produce system- 
atic drifts of the observatories. 

This system will be referred to as CIO terrestrial system, or simply 
"terrestrial system". In the geophysical literature the term "geographic 
system" is customarily used. 

A second two-dimensional frame often mentioned and related to the 
CIO is the following: The axes are on a plane normal to the Xg terrestrial 
axis at the point of intersection with the earth. They are parallel to the 
Xi and Xa axes, but the following sign convention applies; 

X = Xp = Xx 

Y = yp = -Xa 

This frame will be termed "CIO polar system." 

Once the basic reference system has been agreed upon, the mathe- 
matical theory describing the orientation of the instantaneous rotation vector 
of the earth in the "crust-fixed" frame, assuming that body and surface 
torques are zero, is summarized by the followir® matrix differential equa- 
tion {see equation D. 3-13). 

[0} = [Iol[to3+ [iol [to][Iol{a)}+ [w][h}+ {h} (2.4-1) 

where 

[wl = skew-symmetric matrix of the earth rotation vector 
[ lo ] = central earth tensor of inertia of second order 
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{h} s relative angular momentum vector 


Equation (2.4-1) can be solved for the general case of a deformable earth 
by taking into consideration all. the possible variations of [ lo ] and [h] . 
For example, one may write 


[lo] = [I]e+ [AI]r+ [AI]t+ [AI]a+ [AI1p+ other effects (2.4-2) 


where 

[Ih 

[AIIr 

[ATIt 

[AI]f 

[AI]a 

[Al]p 


H initial value of the earth tensor of inertia (independent of 
time) 

s contribution to [ I ] e due to rotational deformation 

s contribution to [ I ]e due to tidal deformation 

H contribution to [I]e due to major feulting (seismic effects) 

= contribution to [ 1]^ due to atmospheric air mass shifts 

= contribution to [ I ]e due to crustal mass displacements 
or plate tectonics 


By the same logic the different contributions to the initial relative angular 
momentum of the earth can be given as 


{h} = [hlg + {h}^ + {h},- + {h]p + {h]p + other effects (2.4-3) 


The values of [Iq] and [h] thus depend on the particular geophysical models 
adopted for describing the earth and its disturbing agents. 

Si effect, one may express 


[AI]p 



[Al]pi 


(2.4-4) 


where n is the number of tectonic plates constituting the earth crust. 
Likewise, 
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(2.4-5) 


> 


[MIt = [aiItj 

•where j = 1,2. . . k and refers to the moon, sun and planets producing tidal 
deformation. 

When all the contributions to the initial values [ I ]e and {h}^ are known 
at some specified times, the differential equations (2.4-1) may be solved 
numerically. To this date an analytical solution of such a complex equa- 
tion seems an insurmountable task, although it is conceivable that in the 
future it may be done by digital computers through symbolic math-logic pro- 
grams. Nevertheless, numerical solutions of first-order differential equa- 
tions are possible by using currently available tested subroutines, such as 
the so called DVDQ -available step, variable order Adanis integrator [ICi'ogh, 
1969], 

Once the equations of motion have been integrated and the vector 
{co} obtained at different time intervals, the coordinates of the pole at these 
specific epochs may be computed by means of equations (2,3-12). 

2.4. 1 Inertia Tensor of a Deformable Earth 

The contributions to the earth inertia tensor due to deformations 
(tidal, rotational, etc.) are obtained through relations involving the external 
potential of the deformed earth and the disturbing potential. 

For example, explicit values for [AT It were derived in [Groves, 

1971] following the above mentioned standard approach. The final equations 
involve the fluid Love number kf and the mass and positions (ephemeris) 
of the external bodies exerting torques. 

By the same procedure Munk and MacDonald [1900, p. 25] give 
the corresponding expressions for [ AT ]r . Nevertheless, recently llochester 
and Smylie [1974] have questioned these results for the first time, based 
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on the tact that the trace of the inertia tensor is not invariant fox’ deformations 
with spheroidal component of degree zero, as was assumed by Munk and Mac- 
Donald, (This is not the case for tidal deformations produced by an external 
body force derived from a potential that is a solid harmonic of degree 2). 

According to [Rochester and Smylie, 1974] equation (5.2,4) in [Mxmk 
and MacDonald, 1960] should be rewritten as (with the notation adopted here) 


[M]r = ^ j^{a)}{w}T-|^tr({to}{co}^-[ArflADl))[l]J (2.4-6) 

where to the first order of {6u} 

tr [ A IfiAo] = 4 Tpe (r ) { x}^ {du} dv (2. 4-7) 

j •'V 

and 

p E (r) s density of the undeformed earth (function of the distance 
to the CM) 

{x3 = position vector of any point P 

[6u3 = vector displacement at point P caused by the deformation 

ks - Love’s number of degree 2 describing the elastic yielding 

of the earth to the second degree term in the centrifugal 
force potential 

G = gravitational constant 

a^ = earth semi-major axis 

Similar expressions for [ AI ]f and [ Al ]r may be consulted in [Sdh- 
chez, 1974] where a different approach through integration is also given. 

Several models for the computation of the contribution to the moments 
and products of inertia due to earthquake faulting or slip displacements 
are available in literature [Mansinha and Smylie, 1967], [Dahlen, 1973], 
[Israel et al. , 1973], but up to now only analytical solutions of the Lagrange- 
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Liouville equations using the approximations established in [Munlc and Mac- 
Donald, 1960, p. 38] have been investigated. 

Thus, the fact that the mass distribution of the earth is subject to 
variations which are functions of time, produces changes in the earth iner- 
tia tensor as expressed by ecpiation (2.4-2). The nature of these changes 
will depend on the earth model adopted and its variations with time. Clearly 
periodic changes in [ I jg may be caused by the effect of [ ATl^ and [ AI];^ ; 
sudden variations, if any, will be due to the effect of [ AI ]f . 

Results reported in [Sdnehez, 1974] show that the effect of the luni- 
solar deformation on io has a maximum periodic variation in the wobble of 
the order of 24 mm. The effects of the rotational deformations are impor- 
tant only in the lengthening of the Kulerian period. 

A widely debated point today is the question of possible excitation of 
the earth rotation vector by earthquakes. This originated from the hypothe- 
sis advanced by Mansinha and Smylie [1967] about the cumulative effects 
of large earthquakes and their influence on the random excitation of the 
wobble. The resulting controversy has continued up to the present time with 
no definitive results established as yet [Mueller, 1975]. Large-scale crustal 
motion that occurs with major seismic activity may abruptly change the 
earth inertia tensor associated to the reference system, and therefore con- 
tribute to a sudden change in the earth rotation vector. 

Schematically Fig. 2. 6 shows the effect of lunisolar tides and sudden 
mass shifts (earthquakes) on the wobble and mean pole, compared to a regu- 
lar curve for the case of rigid earth. 

The earth atmosphere can produce periodic contributions [ AI 
to the initial tensor of inertia due to seasonal variations in the distribution 
of air masses. This can be estimated from observations of pressure on the 
ground. Munk and llassan [1961] showed that the twelfth-month wobble la 
excited atnioapherically, most of the contribution coming from Asia, due to 


the high, pressure over Siberia in. winter. A more detailed discussion in 
[Wilson and Haubrich, 1976] leads to essentially the same conclusion. These 
authors also indicate that contrary to previous beliefs , the motions of air and 
water play an important role in maintainii^ the Chandler wobble. Winds 
acting along the earth surface cause a change [h}A in the earth angular mo- 
mentum. The semiannual and annual components in the spin rotation axis 
of the earth are caused by winds, according to Munk and MacDonald [1960, 
p. 131]. The important finding by Lambeck and Cazenave [1973] was to 
correlate the biennial component {in the rate of rotation) between atmospheric 
zonal winds and the observed astronomical values originally reported 
[Ijima and Okazaki, 1966]. However, today there is still no unanimous 
agreement among the authors as to whether the biennial zonal winds are an 
intermittent phenomenon or their period is 24 months, rather than 26 or 
28 months. Hence a final conclusion cannot be drawn. 

The wobble of the instantaneous rotation axis of a deformable earth, 
after being affected by all the random geophysical effects mentioned above, 
departs from a perfedt conic and is irregular and complex. Other related 
phenomena, such as its excitation, damping or even its period, have not yet 
been fully studied and clarified. 

For example, the rigid earth Eulerian period of about 305 sidereal 
days is extended to what is known, after its discoverer, as the Chandlerian 
period, which is of about 14 months [Chandler, 1892]. As early as 1892, 
Newcomb showed that the lengthening of the Eulerian period is mainly due 
to the yielding of the earth and oceans. 

Chandler [1892] also discovered the twelve-month period, which is 
due primarily to metereological effects, as was mentioned. The composi- 
tion of these two periods (Chandler and annual) produces a beat effect of 
about a six-year period. 

Recently the question of the influence of the liquid core on a possible 
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diurnal wobble has been reopened and theoretically analyzed. The reader 
is advised to consult the papers [McClure, 1973], [Rochester et al. , 

1974] or [Toomre, 1974] where a complete list of classical works on the 
subject is given. 

Fig. 2.7 presents the polhode of the real earth as observed by the 
IPMS during the period 1968.0-1974. 1. Also shown are the barycenters 
(computed centers of the wobble " poles of figure). The determination 
of the barycenters is somewhat arbitrary, due to the complexily of the observed 
wobble. These points will probably not coincide exactly with the axis of maxi- 
mum moment of inertia of the earth at the iiarticular epoch, although they 
should be very close to it. Note that a manifest secular motion of the bary- 
center is present in the figure. 

The secular motion of the mean pole ( = barycenter), if any, will be 
caused by secular contributions to the tensor of Inertia. As mentioned in 
the introduction of this work, one plausible cause, postulated by many au- 
thors, is the effect on [I]f of differential mass displacements due to the 
motion of different tectonic plates in which the earth crust is broken up. 

This being the subject of this investigation, it will be studied at length in the 
following sections and chapters. The next section presents the necessaiy 
formulation for computing the values of the contributions, as used in this 
study. 

It should be emphasized here that the polliodes of the real earth as 
observed by the stations of the international organizations, IPMS (with the 
ILS stations) and BIH, are not coincident, and systematic differences as larg’e 
as 0'^03 are present. 
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Fig, 2.7. Polar Orbit (from IPMS Latitude Observations) for 1962.0- 

1974. 1, Near the Center of the Orbit the Barycenters at Epochs 
1965.0. . . 1971. 0 Are Shown. From [Yumi, 1975] . 








2-5 Total Contribution T All to the Tensor of Inertia Due to Mass Displace - 
ments 


2.5.1 Differential ChanR'es in [ I ] Due to Infinitesimal Changes 
in the Cartesian Coordinates 


By definition the value of the inertia tensor [ I ] may be written in one 
of the following forms 


[I] = 


^[xKxlMm = f px[x][x]^dv 

J n Jv 


where the element of volume is given by 


(2.5-1) 


dv = dxi dxs dxa 


(2.5-2) 


and the density Px is assumed constant for each mass element but not neces- 
sarily equal to the density in the neighboring elements. 

Explicitly equation (2.5-1) may be written 



2 11 

Xa '• X 3 


-X 1 X 2 

Xi+x5 


-XiXa 

-XaXs 

a 

Xi+ Xg 


dv 


(2.5-3) 


The differential changes in [ I ] as a consequence of differential varia- 
tions 6 X 1 (i = 1,2,3) in the mass element Cartesian coordinates are (second 
and higher order partials are neglected) 


[Ml 


where 


d[l] 

Sxi 


i; 6 x 1 =[Mi16xi I- [MslSxg I [At^KjXa (2.5-4) 

1=1 0X1 
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and similar equations hold for the other coefficients. 
Thus, finally 


[AI] = 


f 


0 

t 

-X 3 

1 Px 


2 xi 

0 

/ 

V 


S 

2 xi 


6xi + 


2xji 


-Xi 

0 


0 

~X3 

2xg 


6 X 3 


2 x 3 


0 

2 x 3 


-Xi 

-X2 

0 


6xs 


dv 


(2.5-5) 


2. 5. 2 Differential Changes in [ I ] Due to Infinitesimal Changes in 
the Spherical Coordinates 


The earth being approximately a sphere, it will be more practical to 
obtain the differential changes in [ 1 1 as a function of variations in latitude 
f X , colatitude 60 and radiuS 6r , and at the same time to integrate in func- 
tion of spherical coordinates. 

The transformation between spherical and Cartesian coordinates is 
given by the well-lmown matrix equation 


1 


r sin0 cos X 

X ' 

Xs 1 


r sin0 sin X 

X3 ] 

1 1 

^ rcos0 


(2.5-6) 


and the following linearized differential mapping applies between infinitesimal 
changes in spherical and Cartesian coordinates [Soler, 1976] 


(6X , 60, 6r )- 


-> (6xx,6xa,6xa) 


(2.5-7) 


where the elements of the Jacobian matrix J are 


44 


dXi _ 

ax 

axa _ 

ax 


rsin0sinX= -Xs 
r sin0eos X = xx 


~= reos0cosX= X3C0S X 

otf 


Bxg _ 

a0 


axg _ 

ax 


= 0 


9x3 

B0 


r cos0sin X = XsSin X 

Xi 


= -r sin 0 = 


cos X 


Thus 


X3 


6x1 = “X2 6 X X3COS X 60 I- — 6r 


6xa = Xx 6X+ X3 sinX 60 i- ^ 6r 


6xa = 


_ 3 l _ 50 , ,^2 6 r 
cos X f 


and substituting the above equations in (2.5-5) , follows 


[Ml 


f\ 


2xxXa 

x|-J^ 

X2X3 

1 



-2xxXa 

“X3.X3 

1 


S 


0 

V » 

s 




— 


6X 


2x2XaSinX - 


2xxX3 

cosX 


ar 

ax3 

ax 


= sin0GOS X = — 
r 


= sin0sinX= — 
r 


= COS 0 = — 
r 


(2.5-8) 

(2. 5-9a) 
(2. 0-9b) 
(2. 5-9c) 


-XsXaCOsX-XxXaSinX -x'sCOsX- 


cosX 


„ . 2xx X3 

2XXX3COS X ^ 

^ ^ cos X 


X2X3 

cosX 


->^ainX- 
2xx Xa cos X I- 2x3 Xa s inX 


60 


2 

+ - 
r 


xl+ x| 


-XxXg 

-XxXa 


“X3X3 


xf+x| 


>1 


6r 




dv 


(2.5-10) 


I‘Yom the law of conservation of elements of mass, it follows 


dm = Pjtdv = p |j|dXd0dr (2.5-11) 

(Xf r ) 

where | J | is the Jacobian of the transformation between spherical and Car- 
tesian coordinates. It is known [Soler, 1976] that 

(J| = iHj = hihgha =r®sin6 

Thus, replacing the values of (2.5-6) in (2.5-10) with the above con- 
siderations, finally one has an expression of the form 


[Ml = [M^]6 X+ [Mq] 60 4- [Al,]6r 


(2.5-12) 


where 



a,0,r) 


S0S2X 


-s0c2X 
-s0s2 X 


c0sX 

“C0cX 

0 


r^sin®0dXd0dr 


(2.5-13) 



-s20c^X -S0C0S2X -e20cX 
-s20s^X -c20sX 
o s20 


r^sin6dXd0dr 


(2.5-14) 





1-s® 0 c® X -s® 0s Xc X 
l-s®0s®X 


-s8c0cX 

-s0c6sX 

s=0 


r®sin0dXd0dr 


(2. 5-15) 


In the above matrices, due to space limitations, the following notation equiva- 
lences hold; 

s = sin 
c H cos 
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Observe that as expected a rotation in X does not contribute any 
differential change to the original moment of inertia about the X3 axis. 

When the above Integrals are evaluated over the volume of the earth, 
it is common to assume that the variation of density for neighboring "curvi- 
linear mass elements" is only a function of the radius. Hence 



r) 


- p(r) 


2.5.3 Value of l^ace [AT ] 


Equation (2. 5-12) may be written 

[M] = [ATrot1+ [ATrao] (2.5-lG) 

where 

[AIr„] = [AIq]66 (2.5-17) 

is the contribution to the initial tensor of inertia [ I due to differential 
rotations 6X and 60 while 

= [ATrl6r (2.5-18) 

is the contribution due to expansions or contractions, that is, radial variations. 

It is known that tr [ I ] is invJiriant under rigid rotations of the refer- 
ence system (see section B. 3 in Appendix B). 

Thus, if 

[ I ] E = initial earth inertia tensor 

and 

[ I ] = earth inertia tensor after rotations 6 X and 60 

then 

[I] = [IIe + [AIrot] (2.5-19) 

Therefore, the following implications must hold as a consequence of 


the tensorial properties of [ 1 1 



tr[Il=tr[I]E tr[ AIrot 1 = 0 ==> tr [ AI^^l = tr [ Alg] = 0 


The fulfillment of the above conditions is immediately seen from equations 
(2.5-13) and (2,5-14). 

The situation is different with respect to the contribution due to radial 
deformations as expressed by equation (2.5-lS). In this case after considera- 
tion of (2.5-15) and assuming that the density is only a function of r, one has 


tr[AlRAD] =tr[Al,]6r = 4 /p(r)r5rdv 4 0 

J'i 

which agrees with the previously given equation (2.4-7). 


(2.5-21) 


2.5.4 Inertia Tensor [ I ] in Spherical Coordinates 


From equations (2.5-10) and (2.5-15) it immediately follows 


[I] = 


(2.5-22) 


Thus, 



m=|p 




1-s® 0 cos® X 


■s®6sXcX 

l-s®0s®X 


-s®0cecX 
's9c9s X 
s®9 


r^sin0dXd6dr 


(2.5-23) 


2.5.5 Value of Trace [ I ] 


Using equations (2,5-21) and (2.5-22) it is possible to show that 
(only radial density variations are assumed) 


tr[I] = A + B+ C 




p(r)r®dv 


(2.5-24) 


But in the spherical case 



V = “ irr® => dv = 47rr®dr 

O 


Thus 


/' 


To I 


tr];i] = 8 IT / p(r)r^dr = 8ir r“ / p(r)r'*’dr 

/q *4 


(2.5-25) 


2.5.6 Relative Angular Momentum {h} 
By definition 


{h}= f [x]{x}dm= f p 



but 

ft ft 

-X3X3 + X2X3 

[X]{x}= . 

« ft 

X3X1 - X1X3 

1 

-XgXl + X1X2 

1 


-26) 


_ ft _ ft 

-r^sinX0- r sln6cos0cosXX 

— ft ft 

r^cosXO - r sin0cos0 sin XX 
r® sin^ 0 X 


(2.5-27) 


Thus, finally with 

0 = 60 
X H 6X 


th} = 60 [Ah-}+ 6X {Ah, } 

0 A 


(2. 5-28) 


where 



r^sin0 


(\,0, r) 


-sinX 
cos X 
0 


d0dXclr 


(2. 5-29) 


and 
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} 



sin9 




^ -sin0cos0cos X ] 
-sin0cos9 sinX 
sin®0 


d0dXdr (2.5-30) 


Notice that radial deformations (i.e. expansions and contractions) do not 
contribute to {h}; therefore the following implications are established: 


rigid body ==>> {h} = 0 but [h] = 0 ^=i=> rigid bo<fy 
Finally, the equations given below hold: 

Ah. = AL ; Ah. = -AL ; Ah. = Iga = C 

A-x Am Ag Ai 3 Aa 

2.6 Differential Contributions to the Earth Inertia Tensor Due to Plate 
Motions 


As a real application to the earth of the equations given in section 
2.5, one may compute the differential changes in the earth tensor of inertia 
due to differential tectonic plate motion. This may be expressed by equation 
(2.4-4), namely 

[Al]p = i; [AIrotIpi (2.6-1) 

1 = 1 

where n is the number of tectonic plates constituting the earth crust. The 
elements in the summation of the right-hand side in the above equation will 
be given by 

[AIrotIpi = (6X^[Alj^K + 60,[AlgL) (2.6-2) 

where 

m 2 number of sample blocks on plate 

[Al^Jif s differential changes in the tensor of inertia of block k 
due to an infinitesimal motion 6 Xj. 
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[ = differential changes in the tensor of inertia of block k 

due to an infinitesimal motion 60^ 

Observe that in equation (2.6-1) the assumption fir = 0 is implied; thus ver- 
tical crustal fluctuations are neglected. 

Although conclusive evidence exists supporting the argument that in 
addition to horizontal tectonic displacements of the earth's crust (erogenic 
movements), vertical fluctuations (epirogenic movements) also take place 
(the Fennoscandian uplift being the best documented); the assumption fir - 0 
will continue to be used until precise global data regarding vertical continen- 
tal motion is available, and a unifying framework for vertical movements 
is known. The detection of such movements is hampered by the ditficulfy 
in separating rising of sea level (custatic sea level change) with sinking of 
continents, or vice versa. 

The computation of [Al.jh and[ATnlic as expressed by equations 
(2.5-13) and (2.5-14) involves the integration over the mass of every sample 
crustal block k. Hence changes of density in the adopted earth model must 
be taken into consideration. Therefore in order to perform the integration 
in (2.5-13) and (2.5-14), a properly defined crustal model should be estab- 
lished. 

On the other hand, the summation in equation (2.6-2) assumes know- 
ledge of the plate boundaries and the differential motions for every block k 
(see Fig. 2.3). For this reason some kind of plate model and correspond- 
ing absolute plate velocities with respect to the underlying static mantle 
must be adopted. 

In the following chapter, the particular models adopted in this inves- 
tigation will be described. 
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3 TECTONIC MODELS 


This chapter emphasizes the major assumptions inherent in the se- 
lected models used in the numerical integration and the values of the adopted 
constraints and parameters. 

Two major models are required; a crustal model defining the compo- 
sition of the thin upper layer of the earth, and a plate model describing* the 
boundaries of the different tectonic blocks and absolute angular velocities 
with respect to the underlying mantle. 

3.1 Crustal Model 


Knowledge of the earth’s crustal structure has advanced considerably 
in the last half century, primarily due to investigations in the fields of geod- 
esy (gravimetiry) and geophysics (seismology). 

Seismic research shows a marked boundary betv/een the crust and 
the upper mantle known as th.e Mohorovtcic discontinuily (also referred to 
as Moho or M discontinuity), which separates two layers of very distinct 
density and seismic velocity. 

H: is now fairly well-established tliat the Airy-Heiskanen depth of 
isostatic compensation (partially formulated a long time before the seismo- 
logical data was available) agrees well with the depth of the M discontinuity 
in the ocean basins as well as continental plateaus. This reinforces the 
hypothesis that on a continental scale the earth’s crust is at least approxi- 
mately in a state of isostatic equilibrium. Excess of mass above sea level 







is balanced by a deficiency of mass inside the earth, while the lack of mass 
in the oceans is compensated by an excess of mass beneath. 

Although there are regional disagreements which tend to conform 
with the Pratt-Hayford system of isostatic compensation, these are not 
globally significant and therefore will not be considered. 

To model the upper shell of the earth as closely as possible to phys- 
ical reality in a manner which is computationally feasible, the Airy-Heis- 
kanen isostatie ideas and the formulation given by [Heislcanen, 1938] and 
[Heislcanen and Vening-Meinesz, 1958, p. 137] will be adopted. 

Consequently, the following basic assumptions are made: 

(a) isostatic compensation is complete 

(b) the compensation is local (as opposed to regional) 

(c) the density of the earth’s crust is constant everywhere 
( Pc = 2.67 g/cm*^) 

(d) the density of the underlayer upper mantle (tectosphere) is also 
constant everywhere (p„= 3.27 g/cm®) 

(e) the normal depth of compensation is assumed to be T = 30 km 
Fig. 3.1 shows schematically the cases of emerged and submerged 

blocks with the corresponding notation. 

Written below are the two equations which give the values of the roots 
and antiroots for land and water compartments respectively. 

t=^hjl + 

(3.1-1) 


2T ^ (X, -i l)h . (2T t \,h)[2T I- (X„+ l)h] T(T+ X„h) 
^ —a r." 


t^-uh'jl I i (2T-I- ph^)[2T+ (jLU- l)h] T(T + ph’) (p?-l)h^ ^ 




where 


t s length of mountain root for a land block 
t ' = length of ocean antiroot for a water block (t' > 0) 

h = height above sea level 
h' = depth below sea level (h' < 0) 

Tg = earth radius 


Po = 1,027 g/cm^ S 

sea water density 


Ap 

~ Pn “ Pc 

(3.1-3) 

Xe 

- 

Ap 

(3.1-4) 


_ Po “ Po 
Ap 

(3.1-5) 


For the purpose of this investigation and for reasons explained in the 
following section, an ideal boundary 50 km deep (see Fig. 3.2) which con- 
tains all inequalities of the crust v/ill be taken as depth limit of the upper 
shell of the solid earth model. 

It should be understood that the final results are slightly dependent 
on the type of model selected. This work has used the best model obtainable 
at the present time. Naturally, in the future when a more thorough Imowledge 
of the crustal structure becomes available, changes in the model may be 
expected and consequently the results may improve. 
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Mean oceamc ^ OCEAN p„= 1. 027 g/cm® ' 

depth D~5 Km A ^ 

4- t ^ 


Mean continental 
crust thickness 
T=30 Km 


-i- 


Mean oceanic 

crust thickness CEUST 
To ~15 Km 


JL 

n !■ 

i 


p,= 2.67 g/cm^ 




P=50 Km 


TECTOSPHERE 


pc=3.27 g/crn^ 


Fig. 3.2. Computer Generated Crustal Section Along 8= 60 = © = 30 Obtained Using Equations 
(3.1-1) and (3. 1-2) (Vertical Scale Has Been Exaggerated by a Factor of About 360). 
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Mosaic of plates (adapted from "Plate Tectonics” by John Dewey, copyright ©by Scientific American, Inc. All ri^ts reserved). 

Pig. 3.3. Earth's Principal Tectonic Features. 

(Fig. 2.1 in U. S. Program for the Geocfynamics Project, National Academy of Sciences, 1973). 
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3.2 Absolute Plate Velooity Model 


In order to model the earth’s crustal motions, the most recent geo- 
physical theories are taken into consideration. 

The rigid outer shell of the planet {lithosphere) is divided into plates 
which are in continual motion. Depending on the authors, the number, bound- 
aries and velocities of these plates vary , although there is consistent agree- 
ment with respect to the major important ones (macroplates). 

3.2.1 Boundaries 

The boundaries between the different plates have been defined clas- 
sically by maps of intense seismic activity, even though recently a strong 
correlation between 1° x 1° mean free air gravity anomalies and the plate 
perimeters was reported [Wilcox imd Blouse, 1974]. 

Three types of tectonic boundaries are generally recognized [Cox, 

1973]: 

(i) Ridges : where two plates are diverging, permitting the 
upwelling of magma that creates new lithosphere. 

(ii) Trenches : where two plates are converging, with one plate 
moving beneath the other, eventually to be absorbed into the 
mantle, or "destroyed". 

The direction of relative motion of two plates need not be perpendic- 
ular to the ridge or trenches. 

(iii) l^’ansform faults : where two plates are moving hmgential 

to each other. Lithosphere is neither created nor destroyed. 
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The direction of relative motion of the two plates is exactly 
parallel to the fault. 

Fig. 3. 3 shows the major tectonic features on the earth. 

The plates are assumed to be internally rigid but uncoupled from 
each other . At their boundaries two plates may pull apart or one may slip 
beneath the other, but within the plates there is no deformation. 

In this investigation the tectonic plate model described originally in 
[Solomon and Sleep, 1974] and shown in Fig. 3-4 was adopted. The surface 
of the earth is divided into eleven rigid blocks, avoiding the minor contro- 
versial plates (e.g. , Gorda, Fiji, Iran, etc.) The division between North 
American and South American plates is not yet fully accepted, and therefore 
to combine them in a single plate (American) with a unique set of absolute 
velocities is reasonable. In the figure land masses are drawn up to the lim- 
it of the continental shelf ( 2000 m isobath) while subduetion zones are 
shown as heavy dashed lines. 

Among other plate models available in the geophysical literature 
and previously used by investigators are [Minster et a\ , 1974 and Kaula, 
1975]. 


3.2.2 Ab solute Angular V elocities 

For the purpose of studying the movements of the earth's shell and 
to facilitate treatment of sections of oceanic and continental areas as units, 
it is assumed that the crust is joined rigidly to the upper part of the litho- 
sphere (the tectosphere) and moves with it. This is in accordance with the 
new slab theory which supposes the crust to be a mere passenger driven 
by slabs and clearly contradicts the old concept tliat the continents just float 
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on a viscous mantle. 

Consequently in this investigation crust and tectosphere move as a 
mechanical unit. An imaginary boundary at 50 Ion deep containing the whole 
crust is assumed to approjdmate as closely as possible the effect of the dif- 
ferential displacement of the irregularities of the crust due to plate motion. 
Below that boundary the earth is thought to be structured in some fashion 
not producing any variation in the crustal tensor of inertia. 

Several theoretical models for absolute plate velocities have appeared 
recently in the geophysical literature. Small final differences depend mainly 
on the initial assumptions and the total number of plates considered. 

One of the most complete sets of absolute plate velocities was pre- 
sented in [Solomon et al. , 1975] after postulating that no net torque is 
exerted on the lithosphere as a whole. Absolute angular velocities of the 
Pacific plate relative to the underlying mantle are given, after considera- 
tion of several driving mechanisms which can be modeled quantitatively 
and which are capable of affecting final absolute velocity. 

The description of the different models and their corresponding ab- 
solute angular velocity for the Pacific plate is given in Table 3.1. 

Also given in [Solomon et al. , 1975] are the relative angular veloci- 
ties of the other plates, assuming the Pacific plate to be stationary (see 
Table 3.2). Then tie absolute ar^lar velocity of any plate P^ maybe 
computed as follows 

{^^alpj = + {w^}rp (3.2-1) 


where 


62 




t £^r } Pt = relative angular velociiy vector of any plate 

with, respect to an arbitrary reference plate (Pacific); 

(tCii }rp = absolute angular velocity vector of the reference 

plate (Pacific) with respect to the underlying mantle; 

{ }pj absolute angular velocity vector of any plate Pi 
with respect to the underlying mantle. 

Table 3. 2 tabulates the assumed absolute angular velocity vectors 
}p^ for all the plates (other than the Pacific) computed ti'om the values 
{ tOr]pi and {tOa }rp given in (Solomon et al, , 1975]. 

3.3 Changes (6\fc » 6 Bjj) in Each Block Due to Plate Motions 


Once the absolute angular velocity vector [ }pj for each plate is 
known, the changes and 50 ^ in longitude and colatitude for a sample 
block k belonging to the particular plate may be computed. 

Observe that the components of [ 60^ }p^ in Table 3. 2 can be consi- 
dered differentially small and that the given values conform with a right- 
hand rotation rule. 

The required formulation for evaluating the differential changes in 
the curvilinear spherical coordincates at a point, after a differential rota- 
tion is performed, may be obtained following the general theory presented 
in [Soler, 1976]. 

The final equations are in the form 


6e| 

/ r„ sin 9 cos X 

6X > - 

H"^Il[6m]p^ < i%sin0sinX 

6r ) 

1 r^cosG 


0.3 


(3.3-1) 


TaJale 3.1 

Absolute Velocity of the Pacific Plate* 


Model 

0^1 



Ito 

Orientation of [ to} 

Name 

Description 

10 deg/yr 

Latitude 

Longitude 

A3 

Uniform drag coefficient 
beneath all plates 

-1.490 

2.410 

■6.494 

7.08 

-66°. 4 

121°.7 

B3 

Drag beneath continents 
only 

-1. 800 

3.412 

-8.040 

8.92 

-64*4 

117.8 

B4 

Continents have 3 times 
more drag than oceans 

-1.605 

2.876 

■7.142 

7.86 

-65.2 

U9.2 

C3 

Drag opposing horizontal 
translations of slabs, 
oceanic subduction zone 
only 

-1.333 

2.2U 

■6.147 

6.67 

-67.2 

121.1 

C4 

Same but including ARBT 
and HIMT (see Fig, 3.4) 

-1.353 

2.427 

■6.395 

6.97 

1 

-66.5 

U9.1 

D1 

Maximum pull by slabs 
plus plate drag 

-1. 335 

3.558 

-6.570 

7.59 

-60.0 

110.6 

E2 

Drag beneath 8 mid-plate 
hot spots of [Morgan, 1971] 

-0.755 

2.460 

-6.972 

7.43 ' 

-69.7 

107.1 

E3 

Drag beneath 19 hot spots 
of [Morgan, 1971] ! 

-1. 290 

3.091 

-7.538 

8.25 



-66.0 

112.7 


The components of {wj are referred to the geographic system (x). A right hand rule applies to the vector {to}, 
* From [Solomon et al. , 1975]. 








where 


H s "metric matrix" of the differential transformation 

between curvilinear and Cartesian coordinates; 

R = rotation matrix of the transformation between the 

geocentric and "local moving" frames; 

f ■ skew-symmetric matrix of the absolute angular 

velocity vector for the particular plate Pi contain- 
ing the block k. 

Observe that here [6co]p^ {tOalpj 


For the specific curvilinear coordinates and sign convention used 
in this work, one lias 


H 


r„ 0 0 

r„sin0 0 
, 1 


(3.3-2) 


R= Rg{0)R3(X) 


cos 0 cos X 
- sinX 
sin 0COS X 


cos 0sinX 
cos X 
sin 0sinX 


sin 0 
0 

cos 0 


(3. 3-3) 


and for a given tectonic plate 




0 

- dcca 

StUg 


= 

6t0a 

0 

- 6t0i 



- 6tUa 


0 
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Table 3. 2 

Relative and Absolute Angular Velocily of Plates Other Than Pacific 


Plate Name 

AFRICAN AMERICAN ANTARTIC 


ARABIAN 


CARIBBEAN 


2.406 -4.112 8.973 1.364 -3. QU 5.687 0.675 -3.680 9.596 4.709 -3.784 10.343 1.051 -2.794 6.332 



0.916 -1.702 2.479 

-0.126 -1.201 -0,807 


0.606 -0,700 0,933 

-0.436 -0.199 -2.353 


0,801 -1,236 1.831 

-0.241 -0.736 -1,455 


1.073 -1.901 2.826 

0.031 -1.400 -0.460 

1 

1,053 -1.686 2.578 

O.OU -1.184 -0.708 


1.071 -0.554 2.403 

0.029 -0.053 -0.883 


1.651 -1.652 2.001 

0.609 -1.151 -1.286 

! 

1.116 -1.021 1.435 

0.074 -0.520 -1.861 


•0.815 -1.270 3.102 

3.219 “1.374 3,349 

[ 

•1.126 -0.268 1.556 

2.909 -0,372 2.303 

•0,930 -0.804 2.454 

3.104 -0.908 3,201 

■0.658 -1.469 3.449 

3.378-1.573 4.196 

■0.678 -1.263 3.201 

3,356 -1.357 3.948 

■0,660 -0.122 3.026 

3.374-0,226 3.773 

■0.080 -1.220 2.624 

3*354 -1.324 3.371 

•0.615 -0.589 2.058 

3.419 -0.693 2,805 


■0,439 -0.384 

-0.1G2 

■0.749 

0.618 

-1.708 

■0.554 

0,082 

-0.810 

■0.282 -0.583 

0.185 

■0,302 "0.367 

-0.063 

■0,284 

0.764 

-0.238 

0,296 -0,334 

-0. 640 

■0.239 

0.297 

-1.206 


♦From [Solomon et al. , 19761 
















Table 3.2 (Cont'd) 

Relative and Absolute Angular Velocily of Plates Other Than Pacific 




Plate Name 



COCOS 

EURASIAN 

INDjAN 

NAZCA 

PHILIPPINES 


{Wrlp * -4,715 -14.425 13.332 1.313 -3.569 8.267 6.334 -0.221 10,890 0,693 -9.001 13.692 5.836 -1.362 0.806 

! 


C4 ‘ 


From [Solomon et ai*, 1975], 














After these matrices are substituted in equations (3. 3-1), it easily 
follows that 

6 Xj. = 6 t 03 - ScOiCos XfcCot0fe - ficOgsinX^cot 05 J (3.3-4a) 

60H =-6t0iSinXjc+ fiiCgCos Xfc (3.3-4b) 

When the above values are expressed in radians, one may change 
to linear units as foEows: 


6 Xk = r„sin0fc6Xk 

(3.3-5a) 

60^ = r^60fe 

(3. 3-5b) 


If the plate rotations were not differentially small, the matrix [ 
in (3.3-1) should be replaced by 

Ra(W) - [1] (3.3-8) 

where the rotation matrix Ra('^) is given by (see Appendix A) 

110(^15) ~ costo [ 1 ] - (1 - cosoj) [ P ] + sinco[ a] (3.3-7) 


and 


w = rotation angle along the line of direction 
cosines {of} 

[P] = = projection matrix (3.3-8) 

Notice that for small angular rotations, 

cos to Kd 1 
sin to td 6co 
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) 


and (3. 3-7) reduces to 

Ra(6io) = ll]+6co[a] (3.3-9) 

In this case, 

= 1“' V i = l,2,3 (3.3-10) 

where fitiii (i = 1. 2. 3) are the components of the angular velocity along the 
three Cartesian axes. 

Therefore equation (3.3-9) has the form 

Ra{6t^) = [l]+[^] (3.3-11) 

and finally, expression (3.3-6), as expected, reduces to [ ^ 
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4. 


NUMERICAL EXPERIMENTS AND RESULTS 


Based on the theoiy discussed in sections 2.5 and 2.6 , the differen- 
tial contributions to the earth inertia tensor as a consequence of tectonic 
mass displacements can be evaluated. 

The manner in which this is accomplished practically and the inher- 
ent restrictions will be discussed in the following section. 

Afterwards, the computed values of the matrices involved will be giv- 
en. As a preliminary calculation the earth tensor of inertia is obtained. 
Although there is a possibilily of evaluating [ I ] g approximately for a spher- 
ical earth through integration, using equation (2.5-23) and a density model 
(e.g. , [Bullen and Haddon, 1967]), in this study the latest satellite data is 
used. It should be mentioned here that in this investigation, as will be seen 
later, the matrix [ I ] e bas only a relative role and therefore an approxi- 
mate value is sufficient. 

4.1 Integral Evaluation 

J? 

The basic expressions to be evaluated, namely [ I ] e » [ ^ ] 

A- 

and [ Alg] involve integrations over the total volume of the earth's upper 
layers, under the assumption that the density distribution is known. These 
triple integrals are evaluated along the domains of each plate, so that the 
individual plate contributions can be accounted for. 

The integrations involve two parts: 


(1) Integration with respect to r 

(2) Integration with respect to X and 0 

4.1.1 Integration with Respect to r 

As is known, from the characteristics of the crustal model discussed 
in section 3,1, all variables involved in the integration with respect to r 
are assumed to depend only on the height above or below sea level. 'ITie other 
parameters are adopted constants like T, P and the set of densities. 

Thus a limitation on the integration is imposed by the requirement 
of a global set of earth elevations suitable for digital computation. 

At the present a magnetic tape is available from the Defense Mapping 
Agency Aerospace Center (DMAAC), formerly Aeronautical Chart and Infor- 
mation Center (ACIC), Saint Louis, Mo. , which contains 1° X 1° mean ele- 
vations of the earth solid surface arranged in latitude belts [Czarnecki, 1970], 

It is appropriate to mention now that in areas covered by ice, such 
as Antartica and Greenland, the value given in the tape was used as mean 
elevation of the terrain topography. No other alternative is left to the user 
because no ice elevation is specified. Therefore this work does not consider 
changes in density due to ice coverage. After this report was completed 
the author became aware of the availability of data sets containing Greenland 
and Antartic 1“ x 1“ mean ice thicknesses and topographic surfe.ee eleva- 
tions. This data was compiled by Anderson [1976] from maps supplied by 
the Australian National Antartic Re.search Expeditions (ANARE). Although 
no major changes are expected with respect to the basic conclusions published 
here, a better modeling may be accomplished after the changes in density 
due to ice coverage have been considered. 

Consequently, in this investigation, the crustal layer is supposed 
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to be divided into individual 1 X 1 sample blocks which are not homogen- 
eous, having variable heights and a constant depth of 50 km. 

For every one of these blocks an analytical integration with respect 
to the variable r may be performed. The following two cases are possible 
(consult section and figure 3. 1 for the proper interpretation of parameters): 

(a) Emerged block (h > 0) 

In this case one has 


r,-(T+t) 


p(r)r^dr = Po| r^dr+Pc| r^dr+ p^ 

EMERGED •^r,-{T»‘t) *- 


r^-P 


>„+h)"-[r,-(T+t)]^ 


Ir-(T-^t)]^-(r„-P)^ 


(4. 1-1] 


(b) Submerged block (h^ < 0) 

In a way similar to case (a), it follows 
r^ r„+h' 

/ p(r)r^dr = poj r^dr+ pJ r^dr-> p„ 

*^gj0HERGEO •'I'g+h^ r, -(T-t ) 


r.-(T-tO 


r.-P 


-(r,+hT . . (r,+h')^-[r.-(T-t')]® . _ [r,-(T-t')]^-(r -P) 


These integrals are determined for each block k in a systematic 
computational manner. With the location of the block a search in the earth’ 
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mean elevation tape will immediately identify it as emerged or submerged, 
depending on the sign of h. Entering it into the proper equation (4.1-1) or 
(4.1-2) the integration with respect to r is easily evaluated. 

To account for the earth's ellipticify in a reasonable way, the value 
of r, appearing in the above formulas is the geocentric radius of the earth, 
namely 

r^ - a,, j^l - (f + I f°)cos^0 > I f®cos®0 . j (4.1-3) 

where 

a„ = 6378140 m earth's equatorial radius (4.1-4) 

f = 1/298.257 earth's flattening (4.1-5) 

are the current representative estimates recommended by the International 
Association of Geodesy (lAG) in the fii'St resolution of the XVI*’*' Gcnoi'al Assem- 
bly (lAG, 19751- 

The use of the geocentric radius will introduce what may be considered 
a block deflection, i. e. , each block k will be aligned with its central radius 
vector and therefore will not be normal to the earth ellipsoid. Hie error in- 
troduced is negligible considering that the maximum deflection at about 45“ 
latitude for an ellipsoid the sisie of the earth is only 12', clearly insignifi- 
cant in the context of these computations. 

4. 1 . 2 Integration with respect to X and 0 

The integration with respect to X and 0 requires the plate perimeter 
because evaluation must be done over individual tectonic plates. 

It was mentioned in section 3.2 tliat the plate boundary selected was 
the one given in [Solomon and Sleep, 1974]. Tliose investigators kindly 
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Fig. 4- 1 


Actual Plate Boundaxies tJsed in Evaluating the Integrals 
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pi'ovided a string of data, punched on cards and containing latitude-longitude 
points of the boundary lines. This digitized data was not necessarily ordered 
in a continuous streamline; thus it was preprocessed to make it compatible 
with the information contained on the mean elevation tape. 

Firstly, points at exact 1° intervals (the nearest to the plate bound- 
ary lines) were selected. This set of points was still further reduced to 
evaluate the integration by colatitude belts from X = 0° to X = 360°. Only 
a minimum number of points, all of them significant to the identification of 
each plate, were finally chosen. Fig. 4.1 presents the exact number of 
border points used in the integration process at their corresponding longitude- 
colatitude location. 

Finally, the individual plate contribution to the initial tensor of iner- 
tia due to a differential motion of the crust was numerically implemented 
by using equation (2.6-1). The values (6X^, 60^) were computed for the 
center of each block from equations (3.3-1) which are dependent on the abso- 
lute an^lar plate velocity components. 

Consequently each of the individual 1° X 1° solid blocks forming 
the crustal model has independent differential motions (6Xit, 60^) consist- 
ent with their plate absolute velocity. 

4.2 Earth Tensor of Inertia [I]e 


It is well Icnown that relationships exist between the earth moments 
and products of inertia and the satellite-derived potential coefficients. 
These may be derived by equating the expressions of the earth gravitational 
potential, as given by the spherical harmonic expansion, with the develop- 
ment of the same potential, as expressed by formulas of the MacCullagh 


Thus one obtains 


C's,o = 

|(A+B)-C 

(4.2-la) 

C 3,3 “ 

i(B- A) 

(4. 2-lb) 

II 

O 

E 

(4.2-lc) 

Sa,i = 

D 

(4. 2-ld) 

S3.a = 

F 

(4.2-le) 


Note that the potential coefficients as given a’'Ove have the dimen- 
sions of the moments and products of inertia, that is, ML'\ Nevertheless, 
in practice they are usually in unitless dimensions and generally in a fully 
normalized form. 

Hence the following known general transformations apply tHeiskanen 
and Moritz, 1967, p. 60], 


^ ^ nn 1 

1 

^ = M,a,“ ^ 

1 Cnn 

1 


/2(2n+l)(n-m) ! 
J (n+m) ! 


S'„, J 


^ S„J 

1 

F 

^ S„J 


(4.2-2a) 


Vni = 0 ; = M,a,”Aj2n+l (4.2-2b) 


where 


CuBu Snn = conventional spherical harmonic coefficients of 
degree n and order m (unit L"M) 


Cno> 


unnormalized coefficients (unitless) 


C|ia> SftiD 


fully normalized coefficients (imitless) 


M j = mass and mean equatorial radius of the earth 
From (4. 2-la) and (4. 2-lb) one has 


A = 

C 2,0 ~ 2 C 2,3 '** 

C 

(4. 2-3a) 

B = 

C 2,0 ■*■ 2 C 2,3 

c 

(4. 2-3b) 


Thus finally 


[ Ih = C[ 1] -i- 


C 


^ ~9 ('* 

3,0 ^ 3,3 


S 


“2 S 2,3 3,1 

C 2 , 0 * 2 C 2,3 3,1 

0 


(4. 2-4) 


The value of the maximum moment of inertia may be determined under the 
assumption that the earth has rotational symmetry. 

In this case 


C2,o = A-C (4.2-5) 

ICnowing the dynamical ellipticity of the earth H, as expressed by equation 
(2.3-52) it follows that 


C = 



(4.2-G) 


The value a^ needed in equations (4.2-2) was given previously in (4.1-4). 

The mass of the earth is obtained from the current representa- 
tive estimates of GM^ and G as recommended by the International Associa- 
tion of Geodesy [lAG, 1975): 

GM, = 3.986005 



G = 6.672 s'® kg"^ 

resulting in 

= 5974.228118 X 10®°g 

O o 

The set of fuEy normalized potential coefficients used were those of 
the Smithsonian Astrophysical Observatory (SAO), Standard Earth 111 (de- 
noted SE3) as given in [Gaposchlcin, 1974], namely, 

Cg,o = -4.84170 xlO"^ 

Ca,a 2.3799 X lO"’’ 

Ss,3 = -1.3656 xlO"*^ 

The selection of the SE3 coefficients is such that 

(i) CijO “ Ci^i “ = 0 

(ii) Ca,i = Sa,i = 0 => D = E = 0 

Thus the tensor of inertia [ I ]e will refer to a central system, the third 
axis of which is principal. 

After the selected values are substituted in the general equations 
{4.2-2a and b), one lias 


^3,0 


-1082.637 X 10"^ 

Ca,3 

= 

1.5362188 X 10"*^ 

S.3.3 

- 

-0.88149101 X 10"'’ 

Cs,o 


-0.26311809 xl0^*gcm® 

2,2 

= 

0.37,335410 Xl0^°gcm® 

^ 2,2 


-0.21423268 xl0'‘°gcnT' 


and 


7K 


c 


) 


= 8 - 041506418 X lo'^gcm® 


Thus, finally, recalling equation (4.2-4) the current estimate of the earth 
tensor of inertia [ I ]e , as derived from satellite values is given by 


[IIe 


8.015119938 X lO'^ 0 . 4284653 G 0 X 10 ^° 0 

8 . 015269280 X 10 ^ 0 

( 8.041506418 X 10 '“ 

( 4 . 2 - 7 ) 


where all the units are expressed in gcm^. 

As a consequence of the corollary of postulate 1 discussed in section 
2. 3. 1, it will be proper to assume that the computed center of the earth 
wobble (barycenter) practically coincides with the earth axis of figure, to 
which [ I]e refers. 

Thus one may transform the second-rank tensor of inertia [ I ]e to 
another system with the third axis coinciding with the CIO terrestrial axis 
and the two other axes in the terrestrial equatorial plane. 

Making use of the tensor transformation equation (B. 3-7 ) one may 


write 

[11e ^ 


( 4 . 2 - 8 ) 


where now the rotation angles being small, the value of R is given by 

10 xe 

0 1 -ys I (4.2-9) 

-Xe ye 1 


R - Ilg (-Xg ) Ri (-ye ) « 


and (Xbi y a ) n.re the coordinates of the barycenter in the CIO polar system. 
These values may be obtained graphically from Fig. 2.7. 

For the year 1970 one has 


Xg « 0".028 


(4.2-lOa) 
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(4.2-lOb) 





) 


yg sd 0 .246 


After transforming them to radians and substituting in (4. 2-9), equation 
(4.2-8) finally gives the tensor of inertia with respect to the CIO system 
as (units in g em^ ) 

8.015119938 X10'^ 0.4284053600 X 10'‘° 0. 3033007 '384 x 10^ 

[ I]e = I 8.01526928 X 10^ -0.3129732002 X 10^ 

i 8.041506418 x 10^'^ 

(4.2-11) 

If the matrix [ I ]e is now diagonalized, the principal moments of 
inertia of the earth and the direction of the central principal axis with respect 
to the CIO terrestrial system may be determined. 

After diagonalization, using the eigentheory formalism as explained 
in section B . 4 . 2 of Appendix B, one obtains 


A^p = 8.015108518 X lO'^'^gcm'^ 


B = 8 . 015280700 x lO'^ g cm' 


- 8. 041506418 x 10^ g cnV' 


OP 


(4. 2- 12a) 
(4.2-12b) 
(4.2-12C) 


Tlie spherical coordinates of the positive dh'ection of the central 
principal axes are given by 


Xopl 


X 


op 2 




X piL 

=1 

-14“ 55 '25". 2955 

(4. 2-13a) 

<Ppl 

= 

1 

o 

o 

O 

o 

o 

CO 

o 

to 

Xp2 

- 

75“ 4 '34". 7044 

(4. 2-13C) 

CDps 


0“ 0' 0".2305 

(4.2-l3d) 

Xp3 

= 

-83° 30 '23 ".3089 

(4.2-l3e) 



89“ 59 '59 ".7524 

(4,2-l3f) 
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The direction of Xops was expected, because it corresponds to the 
used values of (4.2-10). 

A major result is the direction of the earth's first principal axis 
which, according to the above values, is practically on the conventional 
terrestrial equator, but about 15° west of the meridian of zero longitude. 

Therefore the earth's central momental ellipsoid will have its great- 
est axis along Xopi and the semimajor axis of a triaxial earth lies 14°. 92 
W of the geographic zero meridian. This agrees very well with the finding 
by 3ursa [1970] that the longitude of the meridian in which is situated 
for a best-fitting triaxial ellipsoid is 14°. 8 W. 

llie transformation of coordinates between the central principal 
and the terrestrial systems may be implemented as follows: 

[x„p3 — {x] or {x] = R{x„p] (4.2-14) 

where 

R = R3{-Xpi)HH(-X0)Ri(-yB) (4.2-15) 

with the meaning of the ^mbols XpijXg and as mentioned above. 

4.3 Crustal Tensor of Inertia f 1 1r 

In this section the results concerning the evaluation of the crustal 
tensor of inertia according to the model described in section 3. 1 are reported. 

The numerical integration of the matrix elements in equation (2. 5-23) 
was performed following the formalism of section 4.1. 

Urn results are presented in table 4.1. Tabulated are the individual 
tensors of inertia for each plate, as well as the corrcs{X)nding total tensor 
of inertia for the whole crust up to 50 1cm depth. 
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Table 4. 1 


Tensor of Inertia for Each Individual Plate 
and the Total Crust {50 Km Deep) 


1 P L A T F 

1 


TENSOR OF INERTIA 

I 

1 NAME ! 

NO 1 


C I 1pi 

1 

AFRICAN 

1 

0.15A9375040D+43 

-.57461126660+42 
0. 36599239270+43 

0*14345257360+42 

0.26123597810+42 

0.38679530540+43 

AMERICAN 

2 

0.45669979430+43 

0.99767866270+4? 

0.37791949010+43 

0.49BB294053D+42 

0.47853782150+42 

0.33954385080+43 

ANTARTIC 

3 

0.32103938X30+43 

-.10982328710+42 
0. 28460049200+43 

0.1<'5022056?0+42 

0.16233108000+42 

0.89046701530+42 

ARABIAN 

A 

0.17796915AOD+A.2 

-.lial730137rM42 

0.16467645060+4? 

-.72e?01?fl940+41 

-.747299^»n4CD+41 

0*24321309360+42 

CARIBBEAN 

5 

0.19538428790+42 

0.5040995'506D+41 

0.29814B4370D+A1 

-.13269803220+41 
0.4572^T96"»6D+41 
0. 20043655140+42 

COCOS 

6 

0,19412A7325D+A^2 

-.18035351620441 

0.977671841PD440 

0.37013024870+40 

0.28377056550+41 

0.19242694290+42 

EURASIAN 

7 

0,33064761550+43 

0.21954097450+42 

0.25444190750+43 

-.38175625690+42 

-.96734867070+42 

0.206745B743D+43 

INDIAN 

8 

0.27358797690443 

0.3761912903D442 

0.16266195640+43 

-.55716172140+42 

0-49073793450+42 

0.29589119770+43 

NAZCA 

9 ' 

0.81882257550*42 

-.41073067850+41 

0.12622831060+42 

-.13013076080+41 

-.22448258630+42 

0.74659768010+42 

PACIFIC 

10 

0*28166704410443 

-.92973343030+42 

0.47876315720+43 

0. 17821657030+42 
-.13439654750+42 
0.49502053700+43 

1 PHILIPPINE 

11 

0. 19303333660+42 

0.1459722513D+42 

0.16358535920+42 

0.69531507100+41 

-.70574419230+41 

0-29852041160+42 






CRUST = [IJc 

1 

1 

0.19765127250+44 

-.16562913680+40 

0.19757955640+44 

-.67571557980+38 
-.A5R75921 37D+40 
0.19011629350+44 




NOTE : ALT TENSORS ARE SYMMETRIC AND REFER TO THE TERRESTRIAL SYSTEM. 
THE UNIT OF EACH ELfHENT ISRorn^ 

H2 


O® POOR QtTALITyi 


Table 4. 2 




Principal Moments and Principal Axis Directions for Each Individual Plate 

and the Total Crust (50 Km Deep) 


1 PLATE 

I 

principal Hor-IENTS 

DIRECTint^ OF PRINCIPAL AXES 1 

1 N A H g ! 

NO i 

r-3-3 

LONGITUDE 

LATITUDE 1 

AFRICAN 

I 

0.A066027917D+A3 

99° 

4 

32.791 

49" 

54' 

49^976 



0.362496^0360+43 

10 8 

36 

36.258 

- 39 

41 

2flo337 



0.13862B9C1QD+43 

14 

42 

37.909 

- 4 

43 

27.103 

AMERICAN 

2 

D.547O2A27O90+A3 

35 

8 

38-160 

18 

1? 

47.466 



0.32507112570+43 

14fl 

15 

5B.R13 

50 

0 

46.551 



0.30206773BBO+43 

112 

14 

7.494 

" 34 

a 

52-129 

ANTARTIC 

3 

0.324502569404A3 

- 14 

?6 

13.690 

? 

25 

52.393 



0.203491P70‘>D+43 

75 

4R 

15.293 

5 

39 

59,093 



0.B6700134550+42 

"127 

33 

56-323 

63 

49 

50.492 

1 

{ ARABIAN 

4 

0.290939M 540+42 

- 24 

10 

6.39B 

- 36 

54 

19.474 



0.28801392370+42 

- 66 

25 

53.872 

43 

38 

24.422 

1 


0.69001599620+40 

46 

28 

5.315 

23 

49 

18.626 

CARIBBEAN 

5 

0.21205n77P0+4? 

140 

33 

40.176 

73 

32 

1.2B3 



0.20945649370+42 

17 

54 

51-556 

9 

3 

34 a 405 



0.412601X6490+40 

105 

41 

53.187 

" 13 

37 

56.729 

COCOS 

6 

0.1973532359D+42 

7 

45 

28-528 

56 

41 

40.086 



0.19528O9234D+42 

" 10 

59 

19.090 

- 31 

53 

17.655 



0.3694234556D+4O 

64 

25 

18.234 

“ 8 

37 

3.376 

i EURASIAN 

7 

0.37183405170+43 

36 

7 

56.223 

- 28 

1 

30*893 



0.29045209300+43 

138 

54 

56.A2B 

" 22 

34 

21 -205 



0.12954845190+43 

81 

54 

56.202 

52 

38 

49.527 

INDIAN 

8 

O.34209S623BO+43 

170 

57 

2.616 

53 

9 

30-414 



0.26311?19P20*A3 

34 

50 

A. 264 

28 

22 

9.125 



Q.1261303090D+43 

112 

35 

0.746 

- 21 

27 

Bi431 

NA2CA 

9 i 

0.B213811173D+4? 

1 

3 

0,632 

" 13 

2d 

44*779 


1 

0.81919290390+42 

" 54 

0 

4.822 

67 

17 

38.B38 



0. 51074553110+41 

86 

36 

17.422 

17 

55 

9.484 

PACIFIC 

1 

10 1 

0.52717755520+43 

113 

4 

57.369 

- 31 

Z 

16.730 



0.4840809?04rM43 

107 

4 

19.431 

58 

49 

20.799 


1 

0.24419226270+43 

21 

26 

45.20B 

- 2 

39 

43.845 

PHILIPPINE 

11 

0.33465004510+42 

37 

32 

53.499 

18 

32 

26.590 I 



0.33206136800+42 

" 93 

42 

4.642 

63 

2 

20.187 



0. 84276741860+40 

134 

7 

30.534 

18 

51 

19.971 


^sss 


Lssas= = 

zzzz 


ssscssssss 

ss...«_SS:3rte 

CRUST 


0.19R1201n76D+4A 

- 88 

56 

55.582 

05 

B 

46.274 



0.1976547396D+4A 

- 11 

46 

47.845 

- 1 

4 

49.400 



0,19757219^30+44 

7R 

7 

50,198 

4 

43 

53.340 
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One should notice that as a consequence of the way the integration 
with respect to X and 0 was determined, the values of these tensors refer 
to the CIO terrestrial system. 

Following the eigentlieory treatment of section B. 4. 2 in Appendix 
B, the symmetric matrices in table 4. 1 can be diagonalized. The resulting 
eigenvalues (principal moments of inertia) and the polar coordinates (X,to) 
of the corresponding eigenvectors (i.e. , directions of the principal axes) 
with respect to the terrestrial system are shown in table 4.2. 

Coniparing the results for the crust at the end of table 4. 2 with the 
values for the earth given in (4.2-12) and (4.2-Ui), the following conclusions 
can be ch-awn: 

1) The principal monients of inertia of the earth crustal layer (50 
km depth) represent only 2% of the principal moments of the 
earth; 

2) The axis of maximum moment of inertia of the modeled crust 
is situated about 5° from the CIO pole in the direction of 
X« - 89°. Recall that the earth axis of figure was situated in 
the direction X - 83°. 5 according to the plotted barycenter 
(epoch 1970.0) as determined from astronomical observations 
processed at the IPMS. 

3) The crust has a principal axis of inertia near the equatorial plane 
but at X'^ - 11°. As known from satellite data, the earth has 

a principal axis in the same plane and at X - 15°. 

A major disagreement with previously reported values appears to 
be in the above results. Milankovich in 1941 [see Scheidegger, 1963, p.l79] 
and Munk [1958] as well as other investigators, have found the pole of the 
earth's "continental-ocean" distinbution near Hawaii. 
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Before interpreting the results of the pole’s location, published to 
date, the reader should be aware of the major differences between the assump- 
tions inherent in the various solutions. These are as follows: 

(i) Milankovich used a 20° X 20° grid, while Munk employed a 

10° X 10° grid in the integration between continental boundaries. 

(ii) Their model was much more simplified. Milankovich's work 
was restricted to standard continental coastlines. He weighted 
the areal mass density for continents differently as compared 
with oceans, without rigorous evaluation of the integrals. 

Munk included in his calculations the continental structure up 
to 1000 fathoms depth. He used a mean elevation of 0.9 km for 
land masses after postulating global isostatic compensation. 

(iii) Both investigators assumed a spherical earth. 

After the discrepancies between the different results were noticed, 
an attempt was undertaken to compare them and new evaluations of the inte- 
grals were made under various modeling hypotheses: 

(a) Crust-M: only the standard crust is considered. (This includes 
the crustal masses contained above the Moho discontinuity and 
therefore excludes the layer of tectosphere down to 50 km depth). 

(b) Crust-M-S: as above for a spherical earth. 

(c) Crust-S: this case uses the crustal model selected originally 
in this investigation, but assumes the earth to be spherical. 

The results concerning case (a) are presented in table 4.3 and clearly 
show agreement with previous findings by Milankovich and Munlc [1958]. 

The maximum moment of inertia of the standard crust is about an 


Table 4. 3 


Principal Moments and Principal Axis Directions for Each Individual Plate 
and the Total Crust (Including Only Down to the Moho Discontinuity) 


I 


PLATE 


J PRINCIPAL MOMENTS 


DIRECTION OF PRINCIPAL AXES 


t 


I NAME I NO I 




LONGITUDE 


LATITUDE 


AFRICAN 

1 

0,1884a97529D+43 

104 

0 

57.391 

45 

59 

48 . 644 



0.1709100992D+43 

100 

20 

24.359 

- 43 

55 

17.557 



0.583Q137625D+42 

16 

15 

22*871 

“ 2 

9 

43.920 

AMERICAN 

2 

0.26985754200+43 

31 

36 

10.913 

18 

16 

59.420 



0*16452430040+43 

134 

41 

42.964 

34 

26 

5.250 



0.14448 703 140 +43 

- 81 ‘ 

20 

2.375 

49 

42 

34.877 

antartic 

3 

0.14458117210+43 

- 13 

9 

36.V39 

1 

50 

41*660 



0.12872301150+43 

77 

1 

13.611 

5 

36 

27.777 



0. 34683539690+42 

r 

“121 

18 

19.767 

84 

5 

41.094 

ARABIAN 

4 

0.16445222160+42 

“ 25 

42 

23.995 

- 35 

2 

18.218 



0.1628719672D+42 

" 70 

41 

16.314 

45 

14 

57.970 



0.37674665420+40 

46 

3 

40.570 

24 

2 

44*226 

CARIBBEAN 

5 

0*10357609630+42 

137 

26 

18.a05 

74 

23 

46.013 



0. 10232427110+42 

17 


42.964 

8 

0 

39.266 



0.2064343915D+40 

105 

16 

19.031 

- 13 

18 

3*898 

i cocos 

6 

0.8006675185D+41 

7 

44 

4 6*306 

56 

49 

34*903 



0,79192759750+41 

- 10 

44 

57.628 

- 31 

47 

49.946 



0.15316130950+40 

84 

33 

23.924 

- 8 

29 

1.797 

1 EURASIAN 

7 

0.21020516030+43 

36 

5 

31.056 

- 28 

37 

6.971 



0.169555362 30+43 

- 40 

2 

49.734 

23 

42 

14*783 



0.68187627940+42 

S3 

10 

35.473 

51 

17 

37.752 

INDIAN 

8 

0.I461500196D+43 

174 

29 

12.996 

53 

44 

45.420 



0.1115313725D+43 

38 

20 

50.218 

27 

52 

5.160 



0.5436504889D+42 

116 

29 

20.108 

“ 21 

14 

15.776 

NAZCA 

9 1 

0*31160649440+42 

3 

3 

4.790 

- 19 

43 

6.794 



0.3101563837D+42 

“ 42 

25 

23.558 

62 

55 

37.535 


i 

1 

0. 19820228460+41 

S6 

26 

57.305 

17 

47 

6.586 

PACIFIC 

10 

0*16099666590+43 

113 

1 

27,209 

“ 30 

43 

32.974 



0. 16875822510+43 

111 

35 

2.580 

59 

15 

58.402 



0.86694577590+42 

22 

38 

53.497 

0 

37 

57.464 

PHILIPPINE 

11 

0. 11621438370+42 

37 

45 

58.055 

18 

28 

31*679 



0.11526136010+42 

“ 92 

47 

9-307 

62 

48 

1.642 



0.30365076620+40 

134 

27 

36.319 

19 

13 

57^492 




CRUST-M 

0-9166113773D+43 

9 

20 

11^869 

“ 14 

45 

4/.450 


0.80607003130+43 

105 

38 

22.271 

" 22 

37 

7.216 


0.85599476430+43 

68 

53 

22.897 

62 

31 

34.217 


Note i TMroNiT of the principal hohemts isgom® 

EAST LONGITUDES ARE CONSIDERED POSITIVE 
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axis with polar coordinates, 

j X = 9°20'll" 869 1 ^ ^ -170° 39 '48". 131 

(<p = -14°45'42".450 ( <p 14°45 '42".450 

This corresponds to an area west of the Hawaiian archipelago, close to 
Wahe Island. The minimum moment of inertia is about an axis passing 
through the northern part of Russia near the city of Gorki. 

The magnitude of the principal moments in this case represents 
only 1 % of the total earth* 

Case (b) gives basically the same results, as can be seen by inspect- 
ing table 4.4. Consequently when only the standard earth crust is integrated, 
the effect of the earth's ellipticity may be considered negligible. This im- 
plies that the asymmetric continental masses {dominant in this case) are dis- 
tributed Over the earth in such a way that their principal axes are not affected 
even when the effect of the equatorial bulge is considered. In general terms, 
this means that land masses are situated primarily north and south of the 
present equator, in a sense contradicting the "Polfluchtlo’aft" theory [see 
Scheidegger, 1963, p. 166] which postulates that continents as floating masses 
on the mantle, will tend to drift to the equator. The above results enforce 
the theory that continents are not just floating masses on the mantle, but 
inhomogeneities forming part of lithospheric plates. 

The last part of this analysis concerns case (c) which introduces 
drastic changes with respect to the original results of table 4. 2 where the 
ellipticity of the earth is taken into consideration. 

Using the crust-tecto sphere model as defined in section 3. 1 but under 
the assumption of spherical earth, the maximum principal moment of inertia 
axis is no longer close to the CIO pole, but at the point with coordinates 
X = - 18° and co •= 10° as presented in table 4.5. 

Consequently the bulge of the earth becomes an important factor when 


Table 4.4 


Rcincipal Moments and Principal Axis Directions for Each Individusl Plate 
and the Total Crust (Including only Down to the Moho Discontinuity) 
under the Assumption of Sphericity 


1 PLATE 

1 

PRINCIPAL MOMENTS 

OIRECTTON DF PRINCIPAL AXES 

1 

1 NAME 

NO 1 


LONGITUDE 

LATITUDE 

1 

! AFRICAN 

1 

0. 18797742 50D+A3 
0 . 1704436039D+43 
0.58270283180+42 

104° 

108 

16 

o' 

21 

15 

0.078 

37.317 

2.708 

45° 
“ 44 
“ 2 

53' 

1 

10 

7^^328 

53.890 

50.195 


AMERICAN 

2 

0.2702458323D+43 

31 

36 

16*348 

IB 

14 

40.570 




0.1649227441D+43 

134 

13 

26*209 

33 

32 

4.2B7 




0.14450905&AD+43 

* 82 

2 

21*3(07 

50 

34 

52.306 


ANTARTIC 

3 

0.1454337697D+43 

“ 13 

8 ‘ 

40*634 

1 

50 

27.806 




0.1295372029D+43 

77 

1 

59.835 

5 

35 

10.133 




0.34B2772750D+42 

“121 

19 

19.188 

84 

6 

59*230 


ARABIAN 

4 

1 

0^16409766240+42 

“ 25 

42 

44*710 

“ 35 

1 

31.810 




0*1625194?B5D+42 

- 70 

42 

27.783 

45 

15 

21.969 




0*37601118550+40 

46 

3 

29*033 

24 

3 

10*598 


CARIBBEAN 

5 

0.10318977420+42 

137 

26 

19.809 

74 

23 

40.352 




0.10194255710+42 

17 

10 

46.538 

6 

0 

41.673 




0,20566969700+40 

105 

16 

21*271 

“ 13 

10 

8.965 


cocos 

6 

0.79738554730+41 

7 

44 

35.641 

56 

49 

2.908 




0*78867797050+41 

- 10 

45 

13*421 

“ 31 

48 

20*462 




0*15256559710+40 

84 

33 

18*766 

- 8 

29 

8*561 


EURASIAN 

7 

0.21061006670+43 

36 

7 

45.372 

- 28 

35 

28*608 




0.16996854870+43 

- 40 

5 

45*455 

23 

35 

57.552 




0. 68242104580+42 

83 

4 

42.967 

51 

23 

46.703 


INDIAN 

8 

0.14586054050+43 

174 

29 

20.074 

53 

42 

37.601 




0,11132326520+43 

38 

25 

3.663 

27 

52 

17*946 



i 

0,54258316770+42 

116 

31 

44*531 

- 21 

17 

0.550 


I NAZCA 


0*31069901780+42 

3 

3 

1*574 

- 19 

42 

4*544 


1 


0*309255S7B0D+42 

- 42 

28 

6.508 

62 

55 

51.125 


1 


0. 19764310760+41 

86 

26 

57*369 

IT 

47 

55.591 


) PACIFIC 

10 

0*1807041551D+43 

113 

0 

34-471 

" 30 

21 

9.043 


f 


0*16839521300+43 

111 

35 

54*136 

59 

38 

23.671 


! 


0*e664118653D+42 

22 

38 

56.977 

“ 0 

36 

55.630 


PHILIPPINE 

11 

0.11580155200+42 

37 

46 

10*956 

10 

27 

43.117 




0.11493056850+42 

- 92 

45 

8*279 

62 

48 

15-946 




0.30261050980+40 

134 

27 

40-898 

19 

14 

25.706 


AftMaacttKic B aaic ] 

1 

KKCOCKesSKEUmCESKXBE 

saaKax 

aaaa 

Kasaasmaa 

aaaaaaaaaavaaaaaa 

aa 

CRUST-H-S 


0.91748924810+43 

9 

5 

55*735 

- 14 

5 

jo.eoe 




0.88677462360+43 

104 

41 

38.107 

“ 21 


32.028 




0. 85503293940+43 

67 

55 

9.757 

64 

7 

55.B17 
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Table 4:. 5 


Principal Moments and I^incipal Axis Directions for Each Individual Plate 
and the Total Crust (50 Km Deep) under the Assumption of Sphericity 


! PLATE I PRINCIPAL MOMENTS DIRECTION OF PRINCIPAL ^ES _ I 

T NAME T NO I C'T'J ~ LONGITUDF LATITUDE T 


AFRICAN 

I 

0.405S39281AD+43 

99° 

3 

42*091 

49 

48 

30.093 



0.36156346580+43 

me 

40 

32*084 

- 39 

47 

33.949 



0.13B3570597D+43 

14 

42 

23*668 

“ 4 

45 

0.949 1 

T 

AMERICAN 

2 

0-5477107520D+43 

35 

8 

42.866 

18 

11 

25.594 



0.32533215970+43 

146 

51 

2*959 

48 

22 

44*691 



0.302535196CD+43 

111 

23 

ID. 819 

- 35 

53 

32.61B 

ANTARTIC 

3 

0,32630642420+43 

- 14 

25 

38-867 

2 

25 

35*103 



0. 28520091620+43 

75 

48 

44*493 

5 

36 

28.667 



0.07040275880+42 

-12 7 

36 

30.814 

83 

51 

21.235 

ARABIAN 

4 

0.29030003480+42 

- 24 

10 

31.322 

- 36 

53 

27.587 



0.28738320370+42 

- 68 

27 

7.456 

43 

38 

54.704 



0.60864354580+40 

46 

27 

53.174 

23 

49 

45.119 

CARIBBEAN 

5 

0.21126562250+42 

140 

33 

36-446 

73 

31 

56.705 



0^20B6flOAQR3D+42 

17 

54 

55.175 

9 

3 

35*239 



0*4112725tt500+A0 

105 

41 

55.776 

- 13 

38 

1.595 

COCOS 

6 

0.19654437380+42 

7 

45 

17.981 

56 

41 

8.625 



0p1944797445D+42 

- 10 

59 

3A.703 

- 31 

53 

47.635 



O.3679019167D+40 

84 

25 

13.232 

- 0 

37 

10.070 

EURASIAN 

7 

0.37255B6922D+43 

36 


13*375 

- 27 

59 

50.431 



0. 29120316850+43 

138 

51 

39.554 

- 22 

27 

4.305 



0.12962011110+43 

81 

47 

25.756 

52 

45 

36.565 

INDIAN 

a 

0.34223 147 2 80+43 

170 

57 

16.565 

53 

7 

18*844 



0.26262319030+43 

34 

54 

46*947 

28 

22 

15.535 



0.12509336990+43 

112 

37 

46.832 

- 21 

.30 

5.907 

NAZCA 

9 

G^Bia9B7Q708D+42 

1 

2 

24.065 

- 13 

26 

17*998 



0.S16S119743D+42 

- 54 

6 

27*613 

67 

10 

33*010 



0.5093077916DM1 

86 

36 

18.585 

17 

55 

57.085 

[ PACIEIC 

10 

0.5263ifl2494D+43 

113 

3 

34.421 

- 30 

44 

57.404 

j 

! 

0*48303375910+43 

107 

1 

57.006 

59 

6 

39.560 

\ 


0*24404442860+43 

21 

28 

42.941 

- 2 

39 

18.730 

PHILIPPINE 

11 

0*33367011720+42 

37 

33 

4.049 

18 

31 

43.843 



0*33108717380+42 

- 93 

40 

15*049 

63 

2 

31.840 


! 

0.84040420120+40 

134 

7 

34.805 

18 

51 

46.694 











1 CRliST-S 


C2.S SS-E ClTS XX M M'S 

0*19763371310+44 

- 16 

22 

3 0.683 

10 

39 

37.666 

1 


0*19780016890+44 

-119 

17 

35.909 

45 

10 

38.267 

j 


0.19771115940+44 

^81 

41 

4*670 

42 

51 

36.037 
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NOTE 1 THE UNIT OF THE PRINCIPAL MOMENTS IS gom= 
EAST LONGITUOFS ARE CONSIDERED POSITIVE 


05> POOR 


the masses below the oceanic standard crust are considered. 

Comparing the results pertaining to the criist in tables 4.2 and 4. 5, 
one actually may conclude that at the present the more dense inhomogenei- 
ties of the earth's crustal layer are situated around the equatorial belt. 
Therefore, the role of the continents as sole agents in the phenomenon of 
balancing the earth's axis of figure lacks consistency. A better intuitive 
picture is attained when the whole plate tectonic structure is considered. 

The layer of higher density below the oceanic plates makes all the difference. 
These masses were systemmatically neglected in previous investigations. 

Hence the present ellipsoidal earth seems to be quasi-dynamically 
balanced, as the distribution of plate masses proves, with the lai'ge dense 
ones (primarily oceanic) occupying the areas around the equatorial bulge. 

Thus the assumption of polar wandering on a large scale in recent geologi- 
cal history seems unlilmly, primarily because the mechanism to explain it 
is not yet available. Nevertheless, a possibility of random secular polar 
motion due, for example, to small lithospheric motions and other redistri- 
bution of matter cannot be ruled out completely. 

Further examination of tables 4. 2 to 4. 5 also shows that the loca- 
tion of the principal moments of inertia axes of each individual plate remains 
practically unchanged, corroborating the fact that ellipticity is not as impor- 
tant as mass distribution. 

Small plates hardly change the direction of the principal axes in the 
case of spherical or ellipsoidal earth. The macroplate change (a few degrees 
only) mostly reflects their assymmetry with respect to the equator. 

A final remark should be made here: If the plates as units are as- 
sumed to rotate without any friction or external forces, the dynamical pos- 
sibility of rotation about their axes of figure always exists. Then the indi- 
vidual pole of rotation will correspond to*the individual maximum principal 


t)0 


moment of inertia axes. 

Although the location of the pole of rotation for the absolute angular 
velocity of the Pacific plate according to table 4. 2 agrees very well in longi- 
tude with the values found by Solomon et al. [1975], for the different models 
(see table 3.1), there is a significant unexplained disagreement in latitude 
(about 30“)- 

4.4 Contribution of Differential Plate notations [ Al]p to the Earth 
Tensor of Inertia 


The contribution [ AT]p to the earth tensor of inertia [ I]e as a 
consequence of differential plate rotations was computed according to the 
theory presented in Chapter 2 and the numerical approach described in sec- 
tion 4. 1. 

In synthesis it will involve the double summation 

t t (6X,[ AT- I, + 66k[AlJ,) (4.4-1) 

i = 1 k = 1 A- o 

H number of plates (eleven) 

- number of 1° X 1“ bloclcs in any particular plate P| 
: given by equation (2.5-13) 

: given by equation (2.5-14) 

: given by eqxiation (3.3-1). It depends on the abso- 
lute plate velocity model. 

Since the expression in (4.4-1) is dependent on the plate's absolute velocity 
modal, several choices are available from the set of eight models charac- 
terized in table 3.1 and given for other than the Pacific plate in table 3.2. 


[Al]p = 

where 

n 

m 

[AlgU 
(6X, 60). 
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Table 4. 6 


Contributions to the Differential Rotation Tensor [ A IrotI 3-s a Consequence 
of Absolute Plate R.otations (Model B4) for Each Plate and the Total Crust 


T PLAT 

E 1 

M 0 

DEL B4 

1 

t NAME 

1 NO 1 

ssxrAixi 

and 

1 

AFRICAN 

1 

“,84839519170+30 

-.29626612960+29 

-.10345841220+28 

0. 30897881820+30 
0.10789778100+29 
0.8494297758D+30 



0*47832573950+29 

-.68403767800+30 

-.47832573950+29 

-.34176817160+29 

0.97869599520+30 

0.0 

AMERICAN 

2 

-.3683 8X004 7D-t-30 

0,13407972060+30 

-.48801027300+29 

0.226334S136D+30 

-.82379025920+29 

0.41718203200+30 



0.14136613110+30 

0.16B4759781D+30 

-.14136813110+30 

-.13028557430+30 

-.3579566736D+30 

0.0 

ANTARTIC 

3 

-.62584424990+27 

0.10924152000+26 
-. 19068178190+24 

0.17924564450+29 

-.31287443590+27 

0.62603493170+27 



0.1162857966D+26 

0. 33299890 lOD+27 
-.1162857966D+26 

-.33315110180+27 

-.19066220060+29 

0.0 

ARABIAN 

4 

0.69177947060+30 

0.7163586124D+30 

0.7418110588D+30 

-.22367284360+31 

-.23162001000+31 

-.14335905290+31 


j 

0. 33913096640+31 

0.11842714240+30 

-.33913096640+31 

0.5188486544D+30 

-.50104632150+30 

0.0 

CARIBBEAN 

5 

! 

0.46636534600+28 

-.20200502450+29 

0.87497989000+29 

-.lA75l48393n+30 
0. 63895696740+30 
-.92161642460+29 


1 

i 

0.440 693 20470+30 

-.90396503130+30 

-.44089320470+30 

0. 68526587490+29 
0. 15620609300+29 
0.0 

COCOS 

6 

0.15199079210+28 

*-.21735695870+29 

0.31083493150+30 

-.62394880010+30 

0.B922883524D+31 

-.31235483940+30 


1 

-.10502043830+31 

0.74725924890+31 

0.10502043830+31 

-.13139558780+30 

-.9X6B074535D+28 

0.0 
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Table 4- 6 (Cont’d) 

Contributions to tbe Differential Rotation. Tensor [A Trot] as a Consequence 
of Absolute Plate Rotations (Model B4) for Each Plate and the Total Crust 


i PLATE 1 

HOD 

E L 

BA 


1 

I N A H E t NO 1 

sexintxi • 

mid 

S60[AlflI 




EURASIAN 


0.68283969970+29 1 136A361020+30 -.36A2617405D+30 

0.10913472880+30 0.60623334100+30 

-.25741869070+30 


-.11110622520+31 0.59076227800+30 -.19019006740+30 

0.11110622520+31 -.11427772160+30 

0.0 


! INDIAN 

B 

-•XB60A33120D+31 

0.67714227810+30 

-.24645963350+30 

-.12371360720+31 
0.45020070570+30 
0. 21068927530+31 



0.2A0110045SD+30 

0.28615200980+30 

-.24011004580+30 

0.230ABA9371D+30 

0.63325216030+30 

0.0 

NA2CA 

9 

-^lB777253LAD+30 

0.65484063790+30 

-.226370069BU+31 

-.27523492570+30 
0.95985825510+30 
0. 24714732300+31 



-,391B17A1A20+31 

0.62703893800+31 

0.39101741420+31 

0.47940559080+31 

0.13746734320^31 

0.0 


PACIFIC 


10 


.11092527060+31 -.58133471060+29 -.10001A62B8D+31 

-.304664612AD+28 -.52415446000+29 

0.11122993520+31 


,33973612010+29 0.32323732580+30 -.38205350060+29 

0-33973612010+29 0.7Z90015051D+30 

0.0 


PHILIPPINE 11 


0.23677655990+30 -.27238027400+30 0.34330045590+31 

0.3133376619D+30 -.39501402000+31 

-.55011422170+30 


0.7759657071D+31 -.10905486790+31 -.30993125390+30 

-.77596570710+31 -.26941912880+30 

0.0 




TOTAL 


-.33718360350+31 0.16667120070+31 10979286470+31 

-.94042641020+30 0.51875548450+31 

0.43122632460+31 


0.59077679290+31 0.12551818210+32 0.47773983650+31 

-.59077679290+31 0.24604695620+31 

0.0 




NOTE I ALL TENSORS ARE SYMMETRIC AND REFER TO THE TERRESTRIAL SYSTEM. 
THE UNIT OF EACH ELEMENT IS R cm^ 


OE POOE OUALM 
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Recently Kaiila [1975] used an independent method to compute absolute 
plate velocities by minimizing the translational motion of tectonic plate 
boundaries. Solomon et al. [1975] model B4 (continents have three times 
more drag than oceans) appears to agree most closely with Kaula's published 
values. According to Kaula [1975] this is the most plausible model, being 
well supported by seismologic evidence. 

For the reasons stated above the B4 model was used initially- in this 
investigation. It will be seen later that all the models lead to the same basic 
conclusions. 

Table 4. 6 contains the values S6X[AT..] and L 60 [ A I _ ] 

i{ A k w 

for each plate and the total crust respectively. It can be noticed that, as 
was explained in section 2.5: 

(a) there is no contribution to [ I]e due to a X rotation. 

(b) the traces of the differential contributions (as a consequence 
of rotations) to the earth tensor of inertia are zero. 

Table 4. 7 tabulates the combined rotation effect. Point (b) ment^ioned 
above still holds if one considers the basic matrix property 

tr([M] + [N]) = tr[M]+tr[N] 

Because the final effect of each plate depends on diverse factors , 

mainly 

(i) Type of plate (continental or oceanic) 

(ii) Dimensions 

(iii) Earth location 

(iv) Absolute velocity 

it is extremely difficult to make an "a priori** forecast of results. 

Observe, for example, that the Nazca plate, although small, gives 
a significant contribution to the initial tensor of Inertia [ I ]e , while the 
Antartic plate, primarily due to location, has a minimum effect. 
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Table 4. 7 


Differential Rotation Tensor [AT rot] = [AT]p for Each Plate 
and the Total Crust Due to Absolute Plate Rotations (Model B4) 


IPLATEI MODEL 04 1 

1 NAME 1 NO 1 [Al]p, 1 

AFRICAN 1 

‘-•60056261770+30 -.71366A 2909D+30 0.27460200100 + 30 

-.488671 5607D+29 0.90940577330+30 

0.84942977580+30 



AHERl CAN 2 

-.2270128736D+30 0.30255569B7P+30 Q.9A04fl93922n+29 

-.19016915040+30 -.4403356995D+30 

0.4171020320D+30 

ANTARTIC 3 

-.61421567030+27 0.34392313300+27 0,17591413270+29 

-.11819261440+26 -.19399094490+29 

0.62603493170+27 

ARABIAN 4 

0.40S3089134D+31 0.83476575490+30 -.17176797820+31 

-.26494 9G605D+31 2 01724 64 ?1 D +3 X 

-.14335905290+31 

CARIBBEAN 5 

0. 44555685820+30 -.92416553370+30 -.78980251790+29 

-.35339521570+30 0 .65477757670+30 

-.92161642460+29 

COCOS 6 

-.10486844750+31 0.74500567930+31 -.75534430790^30 

0.13610393150+31 0.09136956500+31 

-.31235403960+30 

EURASIAN 7 j 

-.104277P2B20+31 0.4-771 1 P6606D + 3O -.55465100790+30 

0. 1300196901D+31 0.49195561960+30 

-.25741069870+30 

j 

INDIAN n 

-.16203230740+31 0.96329420000+30 -.10066511350+31 

-.48656967930+30 0.10035320660+31 

0.21060927530+31 

NAICA 9 

1 

-.41059466740+31 0.692523001SD+31 0.45108210620+31 

0.16344734440+31 0 .2334 53 1 607P+31 

0.24714732300+31 

PACIFIC 10 

-.11432263180+31 0.26510385470+30 -.10303516380+31 

0 .30926965090+29 0.67650605910+30 

0.11122993520+31 

PHILIPPINE 11 

0. 79964336310+31 -.13629209530+31 0.31230733050+31 

-.74463194090+31 -.42195594090+31 

-.5501142? 170+30 

. 


TOTAL 


NOTE : 


0.253593109AD+31 0-lA21Ufs30??D+32 0.2f\79A6971«J0+3: 

-.604019^3390+31 n. 76 A 00 244070+ 31 

0.43122632460+31 

e «====-= *= ==-== ==—=-■ 


ALL TENSORS ARE SYHKD1R1C AND REFER TO THE TFRRFSTRIAL SYSTEM. 
THE UNIT OF EACH ELF Ml NT IS R «'«“ 
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Table 4.8 


Differential Rotation Tensor [AIrqx] s [AI]p Caused by 
Absolute Plate Rotations for All Models 


1 MODEL 1 TOTAL CONTRIBUTION TO EARTH TENSOR OF INERTIA T 

T N A H E T NO F ^AI]. T 

A3 1 

0.95024700520+30 0.1649147382D+32 0.3348E J9 12.00+31 

-.51997094730+31 0.7874B48246D+31 

0. 42494624680+31 

B3 2 

0.47261815520+31 0. 11114399700+32 0.220230A196D+31 

-.91615898270+31 0.75749870700+31 

0.44354082750+31 

3 

0. 25359310940+31 0.142185302 20+32 0.28794697190+31 

-.66481943390+31 0.76480244070+31 

0.43122632460+31 

C3 4 

0.10443500380+30 0,17684984270+32 0.37073485430+31 

-.42395103400+31 0 .76422a3564D*31 

0.41350753370+31 

C5 5 

1 

1 1 

0.71376935170+30 0.16799024510+32 0.35509722780+31 

-.48509224700+31 0,74364603370+31 

0.41371531100+31 

01 6 

0. 1207872A29D+31 0 p 157756692 7D +32 0.333539A159D+31 

-•5256298225D+31 0 .601B229126D+31 

0.4046425797D+31 

EZ ‘7 

j 

0. 21006836290+31 0.14907173250+32 0.40295396870+31 

-.57503107760+31 0.55392909790+31 

0.36496271470+31 

E3 e 

1 

0.35003875150+31 0.12845202320+32 0.30651646540+31 

-.75428931810+31 0.640853Q605D+31 

0,40425056660+31 


NOTE : ALL TENSORS ARE SYMMETRIC AND REPER TO THE TERRESTRIAL SYSTEM. 
THE UNIT OF EACH ELEMENT ISRom" 
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Table 4.8 presents the differential tensor [ AIrotI available 

models. It is evident from the tabulated values that all contributions from 
the different models are very similar. By the common sign of all the prod- 
ucts of inertia in the given matrices, one can predict that the displacement 
of the principal axis of the original tensor of inertia will be in the same 
quadrant for all the models. 

It should be mentioned here tlmt the values in table 4.B correspond 
to a one-year span. Thus they will be the contributions to the earth inertia 
tensor [ I ]e as a consequence of the rearrangement of masses in the earth 
..due to the motion of the plates for a period of one year. 

To this date only a linear value for the plate absolute velocily is 
Imown, and this fact should always be kept in mind. Fortunately, the rates 
of rotation are so small {see table 3.2) that for a period of time ^ 100 years, 
extrapolation of velocities under the assumption of linearity may be considered 
satisfactory. 

4.5 Final Mean Pole Displacements 

The position of the mean pole after the lithospheric plate motions are 
considered may be computed by diagonalizing the matrix [ I ] given by 

[I] = [T1e+[AT]p (4.5-1) 

‘ As implied by the evaluation of [ I ]e and [ AT]p , the tensor [ I ] refers to 
the CIO terrestrial system. In the following analysis, this ideal reference 
frame is assumed to be mantle-fixed, and consequently is not going to be 
affected by the transfer of mass involved in the tectonic fluctuations. As 
mentioned before, this may be considered a reasonable hypothesis for short 
time intervals. 

The differential contribution [Al]p to the original earth tensor 
[ I ]e alters the initial position of the earth principal axes, displacing them 
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after the crustal masses have moved. 

Hence a preliminary problem involves the diagonalization of [ I ] . 

: After this, the new position of the disturbed earth axis of figure will be 

! known, and a comparison with the undisturbed one possible (principal axis 

I before any mass displacement is considered). 

The diagonalization process may be expressed by 

R [[?],+ = N'Je (4.5-2) 

Thus the numerical computation of f i1]e is reduced to the determination 
of the eigenvalues and eigenvectors of the matrix [ I ] given in (4. 5-1). 

Before presenting the results, the relationship between the general 
expression (4. 5-2) and previously given formulae [Darwin, 1877] , [Helmert, 
1880, p. 420], [Tisserand, 1891, p. 485] for obtaining the variation in the 
axis of figure as a result of a given earth mass displacement will be estab- 
lished. 

! 

Selecting the principal axis of the earth as reference system, equa- 
tion (4.5-2) reduces to 

! 

I R[rijE+ [AII pJr^ = riljE (4.5-3) 

I and now 

1 [I] = ^IJE+[AI]p (4.5-4) 

I 

i 

Assuming small coimterclockwise rotation angles between the two 

I principal frames, before and after the crustal motions, one has 
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1 

- 60!3 

fiaa 


60C3 

1 


-6 Ota 
Soil 
1 


[6_a]'+ [ 1] 


(4. 5-6) 
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Thus 


si 


= [[6_ari 

Neglecting second-order terms, 

« [6_a]"[I] + [I][6Q:]+[ll 
and the following conditions are established 


(4.5-6) 


(4.5-7) 


From the above system of three equations with three unknowns, a unique 
set of values Stti (i = 1,2,5) can easily be obtained. For example, using 
Cramer’s rule, the rotation angle about the Eirst principal axis is given by 


6oii => 


-IS 3 

Its 

-Iia 

-Il 3 

Ill“ I33 

I® 

-lie 


Is 3 -I: 

33 ■" l 23 

I12 

“ 

“ Il2 

111” I33 


IlB 


1 2a" Ii 


(4.5-9) 


with similar equations for Stts and Sody. 

A further simplification of the equations (4.5-8) is possible if one 
supposes that the products 6 a 1 1 ^ (V i, j = 1, 2, 3 and i ■/ j) are small, and 
aU of the same order of magnitude. Then their differences aS shown in 
equation (4.5-8) may be neglected. 

Thus (4.5-8) simplifies to 




11360^3 

I3?. ^ 

■ 0 

(4.5-8a) 
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:■ 1 
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(4.5-8b) 
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(4.5-Bc) 
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6cii { Is3- Igg) + 

^23 


0 


(4. 5- 10a) 


6od2 (Iii“ I33) + Iia “ 0 

50^3 (I33" 111) ■*' Iia “ ® 

-whicli gives 

6 ai = (4.5-lla) 

^33^^22 

6 a 2 = - ~ (4. 5“ 11b) 

ill" I33 

60:3 = ' (4.5-llc) 

•*•80" -l-U 

If one Writes 

[I] = mE+[Ai]p = 

(4. 5-12) 
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0 
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a -F -E 

B 0 

-h 

b -D 

1 

<A 

a 
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_s 


(4.5-lOb) 

(4.5-lOc) 


Then substituting the proper values of (4.5-12) in equations (4.5-11) one 
finally has 


6 Oil = 

D 

(4. 5-13a) 

(C+c) - (B+b) 

60 i 2 = 

E 

(4. 5-13b) 

(A+a) - (C+c) 

6oia = 

F 

(4.5-130) 

(B+b) - (A+a) 


The above equations are in the form given by Darwin, and have the disadvan - 
tage that In the case of rotational symmetry, Caia is undefined. 
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IiL this investigation the rigorous approach using eigentheory was 
the selected operational process. Nevertheless, as a consequence of the 
smallness of the rotations, the same results are obtainable from equations 
(4.5-13), if the elements of the matrix [ I ]e given in (4.2-13) and [ AI)p 
from table 4.8 are used in (4.5-12). 

In. the present study the effect of tectonic plate motions on the displac^i- 
ment of the mean pole of the whole earth as well as on the crustal pole were 
determined. 

The results are shown in table 4.9. Several conclusions may be drawn 
from the final tabulated results: 

(a) The first and most important one is that lithospheric motions as 
described by recent geophysical theoretical models of absolute 
velocities [Solomon et al. , 1975] do not produce any significant 
changes in the principal pole of inertia over a time period of 
approximately 70 years. 

This contrasts emphatically with previous published results [Liu et 
al. , 1974] claiming that tectonic plate movement may provide a secular drift 
of the mean pole (10% of the observed value). 

In this respect two points should be clarified: 

(i) The above authors divided the short pa^er In several parts. The 
"Analysis" follows literally the work of Darwin [1B77]. Thus 
primarily equations 4.5-13 are applied. 

(ii) The part of the paper termed "Results of Computations" is de- 
ficient in the presentation of what may be considered to be basic 
intermediary results, such as the integrated values of the prod- 
uct of inertia used in Darwin's equations, or fundamental assump- 
tions like the plate boundary model used. It is coincidental that 
the reported pole direction of secular motion agrees within 2“ 
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Table 4.9 


Displacement of the Polar Principal Axis of the Earth and Crust Due to 
Mass Displacements as a Consequence of Tectonic Plate Motions 



Earth's Inertia Pole Displacement 

Crustal Inertia Pole Displacement 

Model 

1 Year 

70 Years 

1 Year 

70 Years 


Xo{mm) 

Yo (mm) 

Xo (mm) 

yo (mm) 

Xc (mm) 3V (mm) 

Xc (mm) 

yc (mm) 

A3 


-0. 019 

0.56 

-1.33 

0.23 

-1.00 

15. 96 

-69.93 

B3 


-0. 018 


-1.26 

0.14 

-1,01 

9.45 

-70.93 

B4 


-0. 018 


-1. 26 

0.19 

-0.99 

13.37 

-69.30 

G3 

0.009 

-0.018 


-1.26 

0.26 

-0.96 

18.48 

-67.20 

C4 

0.008 

-0. 018 

0.56 

-1.26 

0.25 

-0,94 

17.78 

-66.31 

D1 

0.008 

-0. 014 

0.56 

-0.98 

0.24 

-0.78 

16.94 

-54.74 

E2 

0.009 

-0.013 

0.63 

-0.91 

0.35 

-0.73 

24.43 

-50.75 

E3 

0.007 

-0.015 

0.49 

-1.05 

0.24 

-0,85 

16.66 

-59.43 


Note. The values of the above coordinates are referred to a coordinate sj^stem parallel to the CIO polar 
frame centered at the mean pole XotVo and crustal pole Xc,yc respectively. (See below; sketch 
not drawn to scale). 


yp^Y 


CIO 







with the empirical value published by Markovitz [1960] and based 
on ILS data. 

In the author's opinion the above-mentioned paper has a discon- 
tinuily between Darwin's theory and the final statement that 
plate tectonics may account for 10% of the secular drift of the 
pole, supposedly in the right direction. Practically everything 
in between is left to the reader's imagination. 

If other changes of mass distribution in the earth (eustatic sea level 
changes, melting of ice in Greenland, epirogenic motions, etc.) are postu- 
lated as possible causes of secular drift of the pole, their effect should l^e 
simulated in the most rigorous way before any conclusive answer is estab- 
lished. One feels that the differential contribution to the earth inertia tensor 
as a result of plate tectonic motions is greater in importance, mainly because 
of the magnitude of the masses involved and their asymmetric distribution 
on the earth. Internal shift of masses (i.e., convection currents) is always 
a possibility, although there is no viable way today of knov/ing if they really 
exist. More precise satellite dynamical solutions in the future may clarify 
some questions in this area after reliable data for a number of consecutive 
decades is analyzed. 

(b) While the earth axis of Sgure remains practically unaffected by 
plate motions, even during periods of a century or more, the 
pole of the crust moves only about one centimeter per decade 

in the general direction of the earth maximum moment of inertia 
axis. The vector displacement for the motion of the crustal pole 
is about 50 times larger than the one for the mean pole. 

(c) The general direction of the motion of the pole of the crust coin- 
cides with the prediction of Milankovich and Jardetzky [1962] 
and is in accordance with "Milankovich's theorem" [ifcheidegger, 
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Table 4. 10 

Coordinates Used for the ILS-IPMS Observatories 


J STATION 

PLATE 


C 

□ OHO 

r N A 

T F 

S_ ___ 

INSTR, j 

1 NO 

LOCATiaN. 

COOE^ 


LATITUDE 

LONGITUDE 

TYPE 1 

1 

MTZUSAWA 

7 

39 

08 

03.299 

141 

07 

51. 

FZT 

2 

BLAG0V8SMCH8N$K 

7 

50 

19 

9.300 

127 

28 

45. 

m-i8Q 

3 

CMURtL TOGOT 

7 

47 

51 

56.300 

1G7 

03 

QO. 

V2T 

A 

IRKUTSK 

7 

52 


43 .809 

104 

20 

45. 

2TL"-ino 

5 

MTAR 

7 

39 

08 

01.850 

66 

52 

51* 

ZTL-inO 

6 

nNGELUARDT 

7 

55 

50 

19.966 

48 

49 

GO. 

ITL-inO 

7 

GDkKY 

7 

56 

15 

31. DOG 

43 

59 

00. 

VZT 

8 

POLTAVA 

T 

49 

36 

13.123 

34 

33 

00. 

VZT 

9 

PULKOVO 

7 

59 

46 

15.534 

30 

19 

45. 

ZTF-135 

10 

WAUSAV 

7 

52 

D5 

56.021 

21 

02 

15. 

VZT 

11 

BELGRADE 

7 

44 

4 b 

10.277 

20 

30 

45. 

VZT 

12 

bOlinWlEC 

7 

52 

16 

58-110 

17 

04 

30. 

VZT 

13 

PECNY 

7 

49 

54 

56.000 

14 

47 

15. 

VZT 

14 

DRESDEN 

7 

51 

03 

16.D13 

13 

43 

60- 

TI 

15 

LA PLATA 

2 


54 

31-732 

302 

04 

00- 

VZT 

1 16 

MTLANO 

7 

45 

27 

59.267 

9 

11 

30. 

TI 

17 

TURKU 

7 

60 

24 

57.510 

2? 

26 

45, 

VZT 

IS, 

MT, STRQHLO 

8 

-35 


17.500 

149 

on 

15- 

PZT 

Cl) 

KTZUSAHA 

7 

39 

08 

03.216 

141 

07 

51- 

P2T 

19 

TOKYO 

7 

35 

40 

20.626 

139 

32 

30. 

PZT 

20 

MOSCOW 

7 

55 

41 

56.000 

37 

32 

30. 

PZT 

21 

HAMBURG 

7 

53 

35 

50-100 

10 

01 

00. 

PZT 

22 

NEUCHATEL 

7 

46 

59 

51.620 

6 

57 

30. 

PZT 

23 

GREENWICH 

T 

50 

52 

17.812 

0 

20 

15. 

PZT 

24 

PUNTA INDIO 

2 

-35 

20 

40.500 

302 

4? 

45* 

PZT 

25 

CAGIGAL 

2 

10 

30 

23.0CD 

293 

04 

15. 

PIT 

26 

OTTAWA 

2 

45 

23 

59.850 

2B4 

04 

45. 

PZT 

27 

WASHINGTON 

2 

3B 

55 

16.991 

282 

56 

00. 

PZT 

28 

RICHKONO 

2 

25 

36 

47.020 

279 

37 

00. 

PZT 

29 

CALGARY 

7 

50 

52 

22.000 

245 

42 

30- 

PZT 


KIZIJSAWA 

7 

39 

08 

03.560 

141 

08 

00, 

ASTROLABE 

II-KA-*HI I 

7 

31 

11 

25.000 

121 

25 

45. 

• ASTROLABE 


NOVOSIBIRSK 

7 

55 

01 

46.000 

02 

54 

30. 

ASTROLABE 

(8) 

POLTAVA I 1 

7 

49 

36 

13.123 

34 

33 

00. 

ASTRDLAKC 

32 

POTSDAM 

7 

52 

22 

54.624 

13 

04 

00. 

astrolabe 

33 

GRASSE 

7 

43 

44 

56 *000 

6 

55 

30. 

ASTROLAHF 

34 

6ESANC0N 

7 

47 

14 

56.5B0 

5 

59 

IS. 

ASTROLABE 

C35) 

UCCLE I 

7 

50 

47 

51*890 

4 

21 

30. 

ASTRO LABC- 

35 

UCCLE II 

7 

50 

47 

51.890 

4 

21 

30. 

ASTROLABF 

36 

ALGER 

1 

36 

48 

06-677 

3 

02 

no. 

ASTitOtABE 

37 

PARIS 

7 

48 

50 

09.146 

2 

20 

15. 

ASTRO LABF 

3B 

SAN FERNANDO 

7 

36 

27 

44.000 

353 

47 

45. 

ASTROLABE 

39 

SAD PAULD 

2 

-23 

39 

10.000 

313 

22 

45. 

ASTROLABE 

40 

SAN JUAN 

2 

-31 

30 

39.000 

?9l 

22 

45. 

ASTROLABE 


SANTIAGO 

2 

-33 

23 

57.000 

2B9 

27 

00, 

astrolabe 

(27) 

WASHINGTON 1 

2 

30 

55 

16.991 

2P2 

56 

00, 

ASTROLABE 

42 

QUITO 

2 

-00 

12 

57-000 

281 

15 

00. 

ASTROLABE 


RtCHHONCl I 

2 

25 

36 

47.500 

279 

37 

00* 

astrolabe 

(13) 

PFCNY I 

7 

49 

54 

56.000 

14 

47 

15* 

cz 

* 1 

MIZUSAWA 

7 

39 

08 

03*602 

141 

07 

51* 

VZT 

* 5 

KITAfl 

7 

39 

08 

01.850 

66 

52 

51. 

VZT 

*57 

*58 

CARLDFORTE 

GAITHERSBURG 

7 

2 

39 

39 

08 

OB 

oa.941 

13-202 

8 

202 

in 

48 

44. 

03. 

VZT 

VZT 

*59 

UKIAH 

2 

39 

08 

12-096 

236 

47 

25. 

VZT 

*99 

CINCINNATI j 

2 

39 

08 

19.8 

275 

33 

42. 


* ILS aBSERVATORlES 
' REFERS TO TABLE A- 1 




104 



1963, p. 177]. This implies that the inertia of continental plates 
is less than that of oceanic plates. This is a fe-ct established 
before, and does not contradict any isostatic model. Once again 
global tectonics fits the geophysical scenario perfectly, avoiding 
previously reported contradictions in the direction of polar wan- 
dering [see Munk and MacDonald, 1960, p. 277]. 

4. 6 Apparent Motion of the Mean Pole Due to Station Shifts 

It became evident by the results of the previous section that the earth 
axis of figure (i.e. , the polar principal axis) does not suffer any significant 
displacement as a consequence of tectonic plate motions during a 70-year 
span. Thus, there is no true secular motion of the pole in this period of 
time, at least under the hypothesis and models used in this investigation. 

This section studies the possibility of an apparent secular motion In 
the computed barycenter of the observed wobble, produced mainly by the 
drifts of the observatories monitoring polar motion. 

Previous investigations in this area [Arur and Mueller, 1971] remained 
inconclusive, primarily because only the relative angular velocity for the 
plates available at the time and too noisy observations were used. 

Depending on the selection of the set of observatories, two different 
cases v/ill be analyzed: 

(i) Apparent motion of the mean pole induced by the ILS stations; 

(ii) Apparent barycenter displacements in the IPMS. 

Table 4. 10 compiled from [Yumi, 1975] shows the principal charac- 
teristics of the participating stations. 

^stematic "errors” at the different observatories caused by plate 
rotations were determined according to the following procedure: 
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Table 4.11 

Number of Years of Common Activity for the Different ILS Observatories* 


r 

j Stations Co-observing 

Time Span 

1 

» 

Mizus^wa 

1 

Tchardiui** 

Carloforte 

Gaithersburg 

Ukiah 

Cincinnati 

16 years 

! 

» 

Mizushwa 

Tchardjui 

Carloforte 


Ukiah 

Cincinnati 

4 years 

1 

Mizus^wa 

Tchardjui 

Carloforte 


Ukiah 


4 years 

Mizusp.wa 


Carloforte 


Ukiah 


11 years 

1 

Mizushwa 

Kitab 

Carloforte 


Ukiah 


2 years 

Mizusiawa 

Kitab 

Carloforte 

Gaithersburg 

Ukiah 


11 years 

i 

MizuS^-wa 

i 

Kitab 


Gaithersburg 

Ukiah 


3 years 

1V-T- t 

Mizusawa 

Kitab 

Carloforte 

Gaithersburg 

Ukiah 


19 years 


t Data source: [Dejaiffe, 1972], 

t* The coordinates of Tchardjui and Kitab are assumed to be equal. 

\ 

I 









First, differential changes in the Cartesian coordinates at each sta- 
tion i due to the differential rotations of the corresponding plate Pi were 
obtained using the expression 

[dx]j = [^3 pi[x]j (4.6-1) 

where the vector {x}. was given as a function of spherical coordinates by 
equation (2.5-6); the skew-symmetric matrix [i^co 3pi was defined in section 

3.3. 

After the proper plate rotations are applied to each observatory , a 
new set of disturbed station coordinates is found 

{x}j = {x}j + {dx}j (4.6-2) 

These coordinates define a new Cartesian system (the final or rotated one.) 

The problem is then reduced to determining the angles which will 
transform the original system (x), say at epoch To to the new apparent 
system (x ) at epoch Ti . This displacement, in a sense equivalent to the 
one between the barycenters of the wobble, at epochs To and Ti , is a con- 
sequence of the station motions. Hdnce the position of the (x) system ziiay 
be considered the apparent position of the (x) frame induced by the motion 
of the observatories themselves in the interval Ti - Tq. 

The coordinates of the two systems, assuming differential displace- 
ments, are related by the known transformation 

[x]i = RilxJi = [11] Cx}j (4.6-3) 

where i n the usual notations: 
da-L - e 
Stta = if) 

60 ia 
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are the counterclockwise rotations about the initial (CIO) reference frame. 

By a least square solution, the values of {i = 1,2,3) for each 
particular absolute velocity plate model using ILS and IPMS stations were 
obtained. 

The results are plotted in Figs. 4.2 and 4.3 for a 70-year span 
(Ti “ To = 70 y). The figures show the apparent displacement of the mean 
pole caused by the motion of the ILS and IPMS observatories for the eight 
absolute velocity plate models described in Table 3. 1. 

It is known that the activity of the ILS observatories was sporadic 
and that all the stations were not observing on a continuous basis. 

The discontinuity in the mean pole path in Fig. 4. 2 reflects the changes 
in the number of active cooperating stations throughout the 70-year period 
covered by the observations. The values used are shown in Table 4.11 
where the different groups of co-observing stations and the approximate 
number of years of simultaneous operation assumed, are given. Figures 
4.2 and 4. 3 thus reflect the apparent systematic error in the mean pole 
caused by individual station drifts and accummulated during a 70-yeax span. 
The figures also show the apparent displacement of the equatorial xx axis 
^oint of zero longitude). All the models except B3 (drag beneath continents 
only) and E3 (drag beneath 19 hot spots) produce an apparent westward dis- 
placement at the equator. 

The following conclusions may be drawn: 

(a) An apparent displacement of the mean pole is evidently produced 
as a consequence of tectonic plate motions and should be impli- 
cit in the ILS data as well as the IPMS latitude observations. 

(b) The rate and direction of the barycenter drift depend on the to- 
tal number of observing stations and on the absolute velocity 
plate model. The magnitude of the pole displacement represents 


only 10% to 20% of the astronomically observed value. Model C3 
(drag opposing horizontal translations of slabs, oceanic subduc- 
tion zone only) produces the maximum shift while model D1 
(maximum pull by slabs plus plate drag) gives the minimum dis- 
placement. 

There is not a complete agreement in direction between the 
observed mean pole displacement and the one derived using 
absolute plate velocity models. In general models B3 (drag 
beneath continents only) and B4 (continents have three times 
more drag than oceans) give the best correlation with the ob- 
served pole displacement. 

(c) It is clear comparing the ELS and IPMS results that if the num- 
ber of stations is increased, the amount of displacement is re- 
duced and the direction of the drift changes slightly toward the 
90°W longitude. It also appears that greater numbers of observ- 

. ing stations located on different plates tend to average better 
and produce less appai'ent motion of the ])ole. 

(d) While in the case of possible true displacement of the earth axis 
of figure studied previously, all the geophysical models provide 
practically the same answers; when the apparent position of the 
mean pole is obtained after consideration of station drifts, the 
dependence on the model increases. For exam.ple Model E2 
(drag beneath eight mid-plate hot spots) is in complete disagree- 
ment with the others, producing an apparent displacement of the 
mean pole in a different quadrant ('^10°E longitude). 
Consequently, caution should be exercised in the use of absolute 
velocity plate models, because they are theoretical and do not 
necessarily reflect the actual displacements of the plates . Hence 


Fig. 4.2. Apparent Pole and Reference 
rtiS Stations for DiSerent Absolute 
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Fig. 4.3- Apparent Pole and Reference Meridian Displacement from 
IPMS Stations Ibr DilTerent Absolute Plate Velocity Models 
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these models should not be assumed valid and practically imple- 
mented for the correction of presently observing stations. 

Current observational technology (VLBI, Lunar-Laser Range) pro- 
vides the possibility of establishing a geodetic network capable of detecting 
relative tectonic motions of large plates on the earth with better accuracy 
than the optical astronomical observations. These undoubtedly may contribute 
to the fundamental understanding of global tectonics , providing ground truth 
for comparison with theoretical geophysical models. An initial study was 
conducted [Mueller and Schwartz, 1972] in order to decide if existing stations 
and baselines relative to surface plates could provide the precise measure- 
ments of intercontinental distances needed to recover the plate rates of dis- 
placement. 

To this day no major geodetic results in this area liave been reported, 
primarily because of the limited number of facilities available and the poor 
network geometry. 
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5 SUMMARY AND CONCLUSIONS 


There is conclusive proof now (see Table 1.2) that the mean pole 
of the earth has a secular motion of the order of o". 003/yr toward the direc- 
tion of 70° W latitude. 

The fundamental objective of this investigation was directed to answer 
the controversial question: 

Is there a true secular motion of the pole? 

An examination of the geophysical hypotheses available to explain the 
astronomically observed drift of the mean pole suggest changes in the earth's 
inertia tensor as a plausible cause. 

Rapid accumulation of geophysical and geologic evidence during the 
past decade strongly indicates that the earth surhice is mobile. Large plates 
constituting the outer 50 to 100 km of the earth's crust appear to be moving 
with respect to one another at average rates from 1 cm/yr to as much as 
15 cm/yr. This general structural scheme is termed "global plate tectonics". 

On this basis a differential method was applied to determine the con- 
tributions to the earth's inertia tensor caused by infinitesimal plate motions. 
The approach is consistent for short intervals of time (<100 years), simpli- 
fying the earth model assumptions to present day configuration. 

Nevertheless, the final results are subject to the following modelling 
constraints : 

(1) Crustal model 

(2) Plate boundaries 

(3) Absolute velooify model 


The most recent geophysical ideas have been implemented and used 
in the simulation. Among them, the eight theoretical models for absolute 
plate velocities established by Solomon et al. [1975] provide the framewox'k 
of this investigation and should always count in any final interpretation. 

The major conclusion of this quantitative analysis may be summarized 
as follows: 

No evidence of a true secular motion of the earth axis of figure can 
be found from the earth’s changes of inertia due to tectonic plate movements. 

The indisputable astronomic evidence reflects the apparent motion 
of the wobble's barycenter due to drift of the observing stations as a conse- 
quence of lithospheric motions. 

The above hypothesis is not fully supported by the results of section 
4. 6 where it was shown that the apparent fluctuations of the pole caused by 
the motion of the observatories due to tectonic plate rotations may account 
for only from 10% to 20% of the total observed displacement, depending on 
the absolute velociiy model used and the number of stations involved. 

Thus, neglecting other conjectural influences on the earth tensor of 
inertia due to mass redistributions (sea level changes and ice melting, con- 
vection currents, etc.), primarily because a unifying global framework des- 
cribing the phenomenon is not yet available, the following old question should 
be reopened: 

Is there a sliding of the whole lithospheric crust ? 

The idea preceded global plate tectonics by many years [Wegener, 

1924, p. 152]. Several models have been proposed as possible driving mechan- 
isms, although none of them is fxdly accepted. Among others, convection 
currents in the mantle, mantle plumes or hot spots and, recently, tidel drag 
are the most commonly referenced [Bostrom, 1971 and Bostrom et al. , 1974]. 
Within a short time opposition arose to the theory of tidal drag as a driving 


mechaaism producing westward motion of the crust. Jordan. [1974] in a 
brief note showed that tidal torques are far too small to drag the lithosphere 
at any appreciable velocity. 

A novel, conceptually simple driving mechanism appears implicit 
in the results of this investigation. Having a dynamical base, it fits the 
geophysical as well as the astronomical evidence. 

Assume that the whole crustal upper layer, as modeled in this study, 
can glide over the mantle. According to basic gyro(fynamic theory, the tend- 
ency of the thin shell is to attain dynamical stability, that is, the principal 
axes of the crust will constantly tend to pursue the principal axes of the 
whole earth until they interlock [Munk, 1958], [IhgleB, 1957]. Thus the mechan- 
ism postulated in accordance with the results of section 4. 3 will produce a 
westward drift of the lithosphere at the equator and an eastward displacement 
at the poles, as is schematically shown in Fig. 5.1. 

Crustal gliding over the mantle, as theorized above, is not opposed 
by geophysical formalism. Munlc and MacDonald [1960, p. 282] ix>int out 
that a thin outer shell sliding over the interior could result in a displacement 
of a few degrees at most, and that the stress generated is too small to lead 
to failure. The maximum displacement considered here is about 5°, consist- 
ent with Munk and MacDonald [i960]. The same authors [p* 285] have also 
stated that a crustal gliding over the mantle does not appear to constitute a 
reasonable model for polar wandering; however , as mentioned in the intro- 
duction to this work, a clear difference should be established between secu- 
lar motion of the mean pole and polar wandering. The results of this inves- 
tigation leave the enigma of polar wandering as well as that of large-scale 
continental drift basically unanswered. All the conclusions of this study are 
in agreement with the hypothesis of present global tectonics, not implying 
the support of theories of ]X)lar wandering or major continental drift on a 
geologic time scale. 
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Jardetzlqr [1962] considered the possibility of a gliding of the crustal 
shell over the mantle, analyzing the available paleographic evidence which 
supports the hypothesis: transgressions and regressions, the formation 
of guyots , vertical displacements of individual blocks forming a continent, 
and some features in connection with shear pattern distributions, were re- 
viewed as favorable .evidence in defense of a gliding of the crust. 

Of all the above phenomena resulting from an infinitesimal sliding 
of the crustal shell, the forma.tion of shear patterns may provide the needed 
correlation between real tectonic features on the earth and analytically- 
obtained shear patterns. 

. The first theoretical treatment for calculating such patterns under 
the assumption that the earth crust is ideally elastic, is due to Vening Mein- 
esz [1947]. It is known that the final computed positions of the shear patterns 
with regard to the earth depend on the direction of polar path assumed [Sehei- 
degger, 1973, p. 284]. Recently [Liu, 1974] essentially following Vening 
Meinesz' theory, but assuming the modern values of a shift of the poles over 
80° N along the meridian of 75° W obtained a N-S fracture pattern which 
agrees roughly mth the existing boundary system of the major plates. If 
an equatorial westward rotation of the crust about 14° W is assumed, it may 
provide the other E-W boundaries, such as the Pacific-Antartic ridge or 
South East Indian Rise. 

A westward rotation of the crust is also supported by astronomic 
observations. Recently Proverbio and Poma [1976] considered a group of 
16 BIH observatories located on six different plates. After analyzing time- 
scale data from 1962 to 1967, they obtained a final westward drift of the whole 
crust which accounts for the observed deceleration in the earth rotation. 

Astronomically observed latitude also supports the proposed hypothe- 
sis of the sliding of the crust at the poles in the opposite direction of the 







reduced secular motion of the mean pole. It should be mentioned here that 
Cecchini in an infrequently referenced paper [Cecchini, 1950] after analyz- 
ing astronomical data, shows skepticism about a true secular motion of the 
barycenter, su. lorting the almost absolute fixity of the mean pole. Fig. 

5.2 from [Yumi and Wako, 1970] shows the variation of latitude at the five 
ILS stations for an interval of more than 60 years. The straight lines in 
the figure represent the variation of mean latitude of each observatory if 
the secular variation of the mean pole as obtained Irom the five ILS stations 
is considered. Note that for a sliding crust rotating as mentioned above, 
the variation of latitude will have the trend shown in the figure, with decrease 
or increase in station latitude depending on whether they are in the eastern 
or western hemisphere. It should be noted that the latitude of Carloforte, 
being close to the Greenwich meridian, does not vary greatly. 

It would therefore appear that a crustal sliding over the mantle is 
a possibility that should not be overlooked, and further research should be 
devoted to this topic. 

Fig. 5.2 represents schematically the actual situation at the north 
pole, if an assumed easterly slippage of the crust occurs. 1'wo hypotheti- 
cal epochs Tq '^1900.0 and Tx ~ 1970.0 are shown. Ihe axis of figure of 
the whole earth will hardly change during that period of time. (Although 
clianges in the inertia of the earth due to crustal sliding may move the polar 
principal axis a few millimeters, there is no way to detect this motion wdth 
the present instrumentation). The initial CIO point assumed to be a bench 
mark on the crust wiU move with it, thus causing an apparent motion of the 
mean pole or barycenter, precisely in the direction opposite to the crustal 
displacement. 

In reality the actual apparent displacement of the crust will consist 
in the composition of two movements: 
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Fig. 5.3. Sliding of the Crust and Apparent Mean Pole Positions 



(1) The apparent motion of the mean pole originated by the crustal 
slippage. 

(2) The apparent motion of the barycenter of the wobble due to the 
drifting of the stations, as a consequence of plate tectonics. 

In Fig. 5.4 these two apparent motions are illustrated for a 70-year 
span. One arrow represents the apparent direction of the principal pole 
of inertia due to the eastward motion of the crust {its axis of figure pursues 
the earth principal axis in the opposite direction). The other arrow is the 
apparent motion of the IPMS barycenter caused by station drifts (model 
B4 was used nere). Observe that the resulting vector of these two motions 
closely agrees with the actual observed position of the barycenter. 

All of the above evidence indicates the necessity of taking a new look 
at the secular motion of the pole. 

For a short time span ( < 100 years) the principal polar axis of the 
earth may be assumed to be fixed, that is, mantle-fixed, since it will pierce 
the crust at different points, given the apparent impression that it is moving. 

Tims it is recommendable that all observations be referred to this 
quasi mantle fixed axis. 

It is general practice today to apply what is Icnown as polar motion 
correction to reduce the observations to a common CIO system. Unfortu- 
nately, what may occur is that, depending on the epoch of the reduction, 
different ficticious reference axes are obtained. Hence every polar motion 
reduction today may create a new reference system as a consequence of the 
possible movement of the crust and its attached CIO jjoint. This is the pri- 
mary disadvantage of a crust-fixed reiefenee system: astronomically re- 
duced latitudes (or other coordinates) at different epochs do not refer to a 
common system. 

Tlie principal polar axis of inertia of the earth (axis of figure) has 
the advantage of being dynamical in nature and not striotly geometrical. 
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as the CIO axis. 

It W01 be preferible to refer all observations to a unique mantle-fixed 
system. To be able to do that, geodesists and astronomers need the pole 
wobble referred always to the same point (polar axis of figure) , disregard- 
ing the apparent secular motion of the barycenter . hi this way polar correc- 
tion will refer the observatories to the unique principal axis of inertia of 
the earth. At the same time the international organizations in charge of 
latitude observation should give the secular displacement of the mean pole, 
since this would provide a way of knowing, for example, the rate and direc- 
tion of crust slippage, an important factor in geophysical research. 

After about 100 years a new mean pole should be defined to account 
for possible minor true displacement of the earth’s axis of figm*e. 

Proponents of a mantle-fixed reference system are not few in number; 
among them [Mueller, 1975], [Melchior, 1975] postulated its advantages. 
Previously Fedorov et al. [1972] , after stmfying the effects of relative dis- 
placement of the zeniths at the observing stations, also recommended using 
the mean pole of the epoch of observation instead of the current CIO system. 

Finally, it should be mentioned also that the possibility of a libra- 
tional component is compatible with the above-described theory of crustal 
slippage. 

IE the earth axis of figure has a small periodic motion, as proposed 
by Busse [1970] due to the coupling between the inner core and the mantle 
of the earth, this may create an apparent libration, as seen in the astronomi- 
cally-reduced barycenter. 

In all, a carefiil study about these matters should be undertaken be- 
fore new resolutions on the subject of reference systems are adopted by the 
lUGG or lAU. 
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APPENDIX A 


(Referenced in Sections 2. 1. 3 and 3. 3) 

A. 1 Special Types of Linear Transformations 

lia what follows several important particular cases of the general 
linear transformation L in the Euclidean space if will be studied. 

Consider in general the linear transformation given by 

[x] JU — > {x'} ==?> [x'} = [L]{x} (A. 1-1) 

where [L] is known as the matrix of the transformation and satisfies the 
properties specified in section 2. 1.4. 

Explicitly equation (A. 1-1) may be written: 

= AllXi + AiaXg + A13X3 

Xs = AsiSi + A 233 ^ 2 + A 23X3 (A. 1-2) 

X3 = A31X1 + A3 2X2 + A 33X3 

where 

A 1 j E JK. (M. = set of real numbers) 

The following specific examples of transformation matrices will be 
reviewed* 

[P] s projective transformation matrix 

[Al s antisymmetric transformation matrix 

{ Ra = rotational 
So = symmetric 
[ 1] £ identity 
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A. 2 r*ro1ective Transformations 


When [L] is symmetric. 

and of the form 





[P] = = 

a| 

asQis 


s 

at 

where Od^ (i = 1,2, 3)GJffi. 



and the condition 



3 

tr[to:}[oi}^] = S 

= 1 



1 = 1 


(A. 2-1) 


{A. 2-2) 


is fulfilled, then 

[jP] 

{x] {xH 


(A. 2-3) 


is called a projective transformation and [P] its matrix. 
The above implies 


lxp'} = [P]{x} 


ax(SO!jXt) 
cigdiaiXi) 
olq ( S oilXi) 


(SoiixOfa} 


(A. 2-4) 


Geometrically this means that the point [xp'} is the projection of {x3 on the 
line with direction cosines {a}, hi this respect the following equivalences 
are implied 


s COSQ!i {i*= 1,2,3) (A, 2-5) 

The above geometric interpretation may be shown very easily using classi- 
cal vector representations. 

If (see Fig. A. 2) 
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then 


u 


= unit vector with direction cosines {a} 


OP = txp' } - (u. OQ) . u = 

as was already given by equation (A. 2-4). 

Among the important properties of the projective 
matrix are: 


\ 


transformation 


(i) [P] is idempotent [P]“ ~ [P] 
This is immediately proved, 

(P]2 = [P][P] = (La?)[P] = [P] 

(ii) The matrix [P] is singular 

[P] is singular =c> det[P] = 0 



«! 

OiiOia 

Oil bis 


Oil 

Oia 

0i3 

det[P] = 


r> 

^2 

Oia 0!3 

= Oiittatta 

Oil 

Oia 

0i3 


s 


Oil 


Oil 

Oia 

Oia 


and from the above it follows that 
Rank [P] = 1 = tr[P] 

Therefore the transformation in (A. 2-3) is an application from to 
(i.e. , vectors (points) in the three dimensional Euclidean space E® are pro- 
jected onto the line with {a} direction cosines). 

(iii) The following equalities are easily proved 

{x}^[P][x} = (EdiXi)^ 


and 


{x3^|^[P] + [P'l ] [x] = {-Loc^x.f+iLalx^f 


A. 3 Antisymmetric Transformations 


The linear transformation [L] given by (A. 1-1) is called antisym- 
metric and denoted by [A] if 
[A] 


{x} 5> {xa} [x/ } is normal to [x] 


that is. 


but 


[xa'7[x} = [x}'^[xa'} = 0 


(A. 3-1) 
{A. 3-2) 


Cx/} = [A][x} 

thus one has 

[x 3 ’’{xa'} = [x}^[A][x} = aiiXi h ags^f- (aafl+ a32)XgXa 

E>-ia)^3Xx'i (aifl-i agi)xiX2 - 0 

(A. 3-3) 

Hence, in order to fulfill (A. 3-3) it is necessary and sufficient that 

^11 “ “ <133 “ 0 

^2 3 “ "3^3 

^31 ~ ” %a 

a.1 2 — - 

Using the notation, 

^3 3 ” *131 ■ 

finally one may write 

0 ^ 


[A] s [a] 


as 0 -a^ 

-as ai 0 


and clearly 


[af = - [a] 


Note that 


[x/ } ^ [altx} = tiO - Ixl {a} - - [ a]^ [x] 


{A. 3-4) 

(A. 3-5) 
(A.3-«) 
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and 


which implies in the classical vector representation. 


and 


Xa = a a X where a = vector or cross product 



(A. 3-7) 


(A. 3-8) 


A precautionary note should be made here with respect to the sign eon- 

-> , 

vention implicit in the cross product. The rule for the sign of x a is that of 
a right-handed (or clockwise as viewed along Xa ) rotation which carries 
a into X through an angle < 180 and making sin 6 > 0. 

It is customary in geodesy to use the convention that a positive rota- 
tion is one which appears counterclockwise to an observer looking from the 
positive region of an axis toward the origin. Thus 

“[^1 “ - implies positive counterclockwise rotation 

Fig. A. 1 illustrates the above sign convention with a simple example 




Counterclockwise rotation ; [j] = [i} - j Ak 


-> -> 


Clockwise rotation : [k] [j] = - [i] = IcAj 



Fig. A. 1. Sign Conventions for Clockwise and Countorolookwlso Rotations 
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From (A. 3-6) and (A. 3- 7) it is obvious that in the particular case 
{x} = Ca] 

[x]£x} = [xf Cx3 = -[xltx} = [0] (A. 3-9) 

and 

{xV[xV = {x]"[xl = -[x]"[x]" = COF (A.3-10) 


It is immediate to prove that 

det[^] = 0 and Rank[^] = 2 
Notice that if 


then 



a ^ u (unit vector) 



■ 0 

— 


[u] = 

Oia 

0 



--0^.2 


0 __ 


[a] 


where {a} are the direction cosines of the unit vector u. 


A. 4 Orthogonal Transformation s 


It is very well Imown that in general a transformation like (A. 1-1) 
is called orthogonal when the matrix of the transformation fulfills the properly 

[Nf = [N]-^ (A.4-1) 

If det [N] = +1 the transformation is called "proper", or "direct". Every 
orthogonal direct transformation conserves the length and the angles between 
vectors. 


A. 4. 1 Relations between Orthogonal and Antisymmetric Transfor- 
mation Matrices 

(a) Direct orthogonal matrix as a function of antisymmetric matrix 
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It may be proved that if [a] is the matrix of an antisymmetric trans- 
formation, then 

[N] =[[1] -[a]j [[l]+[a] j-" (A. 4-2) 

is the matrix of a direct orthogonal transformation. 

Observe that det[^] = 0 bnt det [[ 1] + = 1 5 ^ 0 

Given a square matrix [^] one- has 

[[a] + 11]] [[a] -[!]]= laf-in - [ [a] - [1]] [[a] + [ 1]] 

In the same way it may be proved that 

[[1] - [a] ] [ [1] + [a] ] -^ = [ [1] + [a] ] -^[[l] - [a]] 

Then, if [a] is antisymmetric, and 

[N] = [[!]-[ a]'] [[ll + [a]]'" (A. 4-3) 

one has 


[Nf =[[!]+ [a]^ -[a]^]=[[l] - La] ]-"[[ll+ [a]] 


but 


thus 


[N]-" =[[!]+ [aj] I’ll] - [alj-’- 


[N]' = [N]-^ 

and the statement related to equation {A. 4-2) is proved. 


(b) Antisymmetric matrix as a function of a direct orthogonal matrix 
Reciprocally, if [N] is direct orthogonal, then 


[a] =[[N]+[l]]'"[m-[N]J 

(A. 4-4) 

and the following condition must be fulfilled 


det [[N]+ [1]]t^ 0 

(A. 4-5) 

Let equation (A. 4-3) be written in the form 
[ [11+ [a] ] = [1] - [a] 

(A. 4-6) 


[N][l]-i- [N][a] = [l]-[a] 
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or 


which gives 


[[Nl + Ul] [a] = [11- [N] 

Consequently, if condition (A. 4-5) is fulfilled 
[al =[[N1 + [1]]*^[[1] -[N]] 

which proves equation (A. 4-4). 

Now it will be proved that the given expression for [^] in (A. 4-4) 
is antisymmetric. After transposing (A. 4-6) 

[[l] + [ar][Nf = [l] -rar==>[l]+[ar= [[1] -[ar][N] 

(A. 4-7) 

Multiplying (A. 4-6) on the left by [N]’^ 

[l] + [al - [N]-"[[l]-[a]] (A.4-8) 

Multiplying {A. 4-7) on the right by [ 1] [a] as given by (A. 4-8), one has 

[[11 + [a]^] [[11+ [al ] = [[1] - [ar][[ll - [al] 
thus finally, 

[al + [a]T [al - [a]" i 

2 [[ar+[a]] = [01 => [ar = -[a] 

A. 4. 2 Orthogonal Transformation Matrix for a Rotation around 
a Line of Direction Cosines [a} 

In the following treatment all counterclockwise rotations as defined in 
section A. 3 are going to be considered positive. 

Suppose one has a vector [x] == OQ (see Fig. A. 2) referred to a 
Cartesian frame (x ). A line of direction cosines [a] is given and a rota- 
tion 0 around this line performed. 

After the rotation the vector [x] will be transformed into the [xr }. 
The problem is to find the orthogonal matiMx Ra(S) of the transformation. 


13 1 



tx} 


Ha{0) 


-> [xO or (xrI = Ra( 0){x} 


(A. 4-9) 


In Fig. A. 2 the vectors OP = [xp} and OA s {x a } transformed from 
{x} by the previously studied projective and antisymmetric transformations 
are shown. 

That is, 

OP = [x'p } = IP]{x} (A. 4-10) 


where, recalling equation (A- 2-1) 
[P] = 


and 


but 


OA - {x/} = [a]^ [x] = [X] [o'} 
From the picture, clearly 
[xh } = or = op + ^ 

ife = OM cos 6 - sin 0 


(A. 4-11) 


(A. 4-12) 
(A. 4-13) 


and 

(?M = 6q - OP (A. 4-14) 

Therefore, replacing the value given by (A. 4- 14) in (A. 4-13) and this result 
in (A. 4- 12), one has 

[xr } = OP + ( - OP) cos 0- OA sin 0 (A. 4-15) 


Substituting above the values from (A. 4-10) and {A. 4-11) one obtains 

[xrO = [P]{x] + I [x] - [P][x] I cos 0+ [a] [x} 
or 

[xr'}= |^[P] + [[1] - [P]] cos0 [a] sinej[x] 

Hence, finally recalling (A. 4-9) 

lia (0 ) - cos 0(1] > (3 - cos 0)( 1’ 1 + sin 0 [ a 1 
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(A.4-l«) 


which gives the value of the rotation matrix R a (6) as a function of the known 
parameters; the rotation 9 and the direction cosines {cd}. 

Note that 

rJ( 0) = Ra(“0) = cose[l] + (l-cos0)[Pl“Sin0[a] (A. 4-17) 

Clearly the above equation will apply to the clockwise rotation case. 

In order to prove that Rj (9) - Ra^(9) (observe that this cannot be 
done by direct inversion of (A. 4-16) because [P] and [a] are singular matrices) 
one will show that 

Ra (0)Ra(0) = 1 
After multiplying the matrices 

W(0) = [P] + [l]Gos0 - [P]Gos0- [alsinO 

Rg (0) = [P] [I]cos0 - [P]cos0 [a]sin0 

and after simplification, considering that 

(i) [P] is idempotent => [Pf = [PI 

(ii) [a][P] = [P][a] = [0] 

one has 

Ra (0)Ra(0) = [PI + COS®0[1] - cos®0[P] - sin®0[a]® 

It may be proved that (see B. 1-9 in Appendix B) 

[al^ = [P] - {tr[P])(l] 
but tr[P] = 1, thus 

[a]^=[P]-[l] (A. 4-18) 

Then finally 

Ra(6)Ra(0) = [P]+ cos®0[l] - cos®0[P] - sin®0[P] + sin®0[l] = [1] 
and the orthogonality condition of Ra(0) is fulfilled. 

Remark . The transformation matrix Ua(0) as derived above rotates {coun- 


terclockwise) vectors in Rp while keeping the basic frame (x) fixed. 


That is, it transforms the vector components (w.r.t. the same fixed Carte- 
sian system) before the rotation to the components after the rotation. In- 
stead, when the reference frame is rotated (countercloclcwise) while keeping 
the vectors fixed, their components will be transformed by the matrix Ra{6). 
The basic concept is self-explained below in Fig, A. 3. 

a) Rotation of Vectors (Frame- Fixed) b) Rotation of Frame (Vector- Fixed) 



Ra(6) : Counterclockwise rotation of vectors <^> Clockwise rotation of frame 
Ra(0) : Clockwise " " " <=> Counterclockwise ” " 

Fig. A. 3. Rotation of Vectors and Frames 

Ih particular, the following notations apply when the counterclockwise rota- 
tions are taken around the three Cartesian axes: IIqi (0) ?? Ri (0) V i = 1, 2, 3 
where o?i(i = 1,2,3) are the direction cosines of the xj Cartesian axes. This 
can be easily proved using (A. 4-17). 

A. 4. 3 Knowing the Orthogonal Matrix Ra(6)j to Obtain the Line 
about which The Rotation 0 Is Performed 

The problem will be divided into two parts: 

(a) Computation of the direction cosines Ca3 

As shown in the previous section, every orthogonal rotation matrix 
keeps the direction cosines (o:] of the rotation axis invariant. 

This Is also true for the matrix Ra(6). 





Thus, for any vector [x} of {a} direction cosines, it follows that 
Eq{0) tx} = Ra^(0){x} = [x} 

whence 

tRa(e)“V(0)] W = [a]{x] = [0] 

which implies 

— > — > “> — > 

a A X = 0 => a and x are collinear 
Thus the components of a will define the axis of rotation, that is, 

0 ria-rai rig -131 

r2i-ria 0 

rai-rx3 rgg-rgg 0 

If the counterclockwise sign rule as defined in section A. 2 is followed, 


[al = Ra(0) - Rj(0) = 


8-1 " ^33 “^as 


83 =r 3 i“ris 
83 ” ri g-rgi 

and the direction cosines sought will be given by 


(A.4-19a) 

(A.4-19b) 

(A.4-19C) 

(A. 4-20) 


(b) Computation of the rotation angle 0 

It is known that the trace of a second-order tensor is invariant under 
rotation (see section B.3. 1). Therefore from equation (A. 4-16) it follows 


and 


tr[Ra(0)] = 3COS01- (l-cos0) = 1+ 2eos0 


COS0 


tr[Ra(0)] - 1 
2 


(A. 4-21) 


As an illustrative example, assume that the following orthogonal 
matrix is given 
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10 0 
0 cos 6 sin0 

0 -sin0 COS0 


(i) Direction cosines of the rotation axis 
From (A. 4-19) it can be found 
ai = 2 sin0 

ag = 0 

aa = 0 


Thus 


Ofi = 


2 sine 


4 sin® 9 


= 1 


ffg = CHa = 0 


(ii) Angle rotated 
Using (A. 4-21) 


_ 2 GOS0 1-1-1 Q 

COS0 = r = COS0 


Therefore, the given matrix is, as expected, Ri (0) : rotation of an angle 
0 around the first axis. 


A. 4. ^mmetric (Orthogonal) Transformations 

lii order that Ea{0) be symmetric, it is necessaiy and sufficient that 

sia0[a] = [0] 
which implies 

sin6 = 0 0 = niT VnGlL 

Two particular cases may be studied: 

(a) n is an even number, n = 2k , VkGJE. 


Then 


0 - 2kir and Ra (Sktr) = [ 1] 
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(A. 4-22) 




) 


This was obviously expected, because wheu a rotation of a multiple of 2 tt 
around a line is given, the transformed vector is itself. Thus the unit matrix 
is a particular case of rotation for 0 = 2 kff . 

(b) n is an odd number n=2k+l ,Vk GJffi. 

In this case, 

Ra(7T) = 2IPl -[1] 


and 


thus 


Ug] = ■Ra{ir){x} ==^ {xs'l = 2[P]{x 3 - [l]{x] 


tx/}+ {;x} = 2[P1Cx3 

which implies that [xs'} is symmetric to {x] with respect to the line of direc- 
tion cosines {o;}. 

Therefore, the matrix of a symmetric transformation [S] is given by 
Se = Ra{7T) = 2[P] - [1] (A.4-23) 

Thus the matrices of the symmetric transformations are orthogonal 
and symmetric. 

Clearly the symmetry with respect to the origin may be obtained 
from (A.4-23) when [a] = {o} . 

-10 0 

-1 0 

s “1 

In particular, symmetries with respect to the three axes are easily 
obtained. For example the symmetry respect to the first axis is given as 

10 0 

-1 0 

s 

In general, denotes the symmetry matrix whose elements s nn 
are equal to zero except 


So = -[1] = 


Si = 


*11 


- I 


1 ~ Sjfk 


-1 
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APPENDIX B 


(Referenced in Sections 2.2, 2.5.3, 4,2, 4.3, A. 4, D. 1, D. 2 and C. 2) 

B. 1 General Comments 

This appendix summarizes the transformations between the second- 
order inertia tensors [ I ] and [ I ] as defined in Section 2. 2. When the con- 
text makes the difference between the symbols [ I ] and [ E ] clear, the word 
"associate" for the inertia tensor will be omitted. 

The accuracy of present observations (VLBI, Lunar Laser Range) 
will require consideration of these tensor transformations between selected 
frames of reference. It is standard practice today to give the principal mo- 
ments of inertia of the earth without clearly specifying the system to which 
they are practically referred. The same situation applies to the set of satel- 
lite-derived potential coefficients. 

This work treats the tensors of inertia as individual entities. This 
generalization results, in particular cases relating moments and products 
of inertia, some of which may be found in classical textbooks on dynamics 
[Routh, 1905, I]. From the point of view of the reader, the anlysis is greatly 
simplified because only the basic concepts in matrix algebra are needed to 
establish the dynamical fosmalism. 

Based on Binet's theorem, a new approach is presented for computing 
the moments of inertia and the directions of the principal axes at any point 
of a body. This is generalized to the diagonalization of matrices using par- 
ticular rotation matrices. 

B. 1. 1 Relations between the Tensors [ I ] and [ I ] 

The following matrix equality Is easy to ]>rove, 
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= [yl[xl+ ({x]Uy])[l] 


(B. 1-1) 



where the coefficient in parenthesis may be expressed in one of the following 
forms: 


{x}'Cy}= ty]Ux] = tr[{x][y}'] = tr[ty}tx}"] (B.1-2) 

Similarly, 

{yltx}" - [[x}{y}^r = [x][y]+ ({x}n.y}) [1] (B.1-3) 

and properties analogous to (B. 1-2) hold. 

From (B. 1-1) and (B, 1-3) one has 

{x} {y}^ = [yllx] - [x][yl + [y} [x] ^ 
or 

[y}{x}’' = rx][y] - [y][x] + [x}{y}^ 

Equations (B. 1-4) or (B. 1-5) immediately yield 

[xH_yr = [yl'[x]+ {xHy}'- [yjW (B.1-6) 

Consequently, in general 

[x][yr 7^ [y]"[x] 

In the particular case when {x} = {y } the following matrix rela- 
tions may be established: 

(i) From (B. 1-6) 

[x][xr = [x]'[x] 

which is obvious from the basic definition of skew-symmetric matrices and 
was already given by equation (2.2-2). 


(B.1-4) 


(B.1-5) 
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(ii) From {B. 1-1) 

{x}{x}' = [x]^+ ({x}"{x})[l] 


{B. 1-7) 


where clearly recalling (B. 1-2) 

{x3^[x] = tr[{x]{x}’^] (B.1-8) 

Thus finally, one concludes 

Ix]® - Wtx]"- (tx}Mx})[l] {B.1-9) 


The above equation gives another way of expressing the associate ten- 
sor of inertia, 


[II 


{ {x3 [x]) [ 1] - {x} [x}^jdm 


(B. 1-10) 


Taking into consideration (2.2-4) and (2.2-5) the relation between the 
matrix tensor [I] and [I] follows immediately, 


[II = (tr[I])[l] - [I] (B.1-11) 

Applying traces on both sides of the above matrix equation, one has 
tr [ I ] = 3 (ti- [ I ] ) - tr [ I ] 

which results in the relation between the traces of [ I ] and [ I ] 

tr[I] 2tr[Iil (B.1-12) 

Hence it is now possible to express [ E ] as a function of [ I ] . Introducing 
(B. 1-12) in (B. 1-11) 

[I] = ^(tr[I])[l]-[I] (B.1-13) 


From (E.1-13) or (B.1-11) it immediately follows, 
A < B < C 111 ^ > I 33 
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B-1.2 Variatioiis of [I] and [I] with respect to Time 


The following notation will be employed 


[f] 


d 

dt 


m 



where as usual the dot means differentiation with respect to time. Similar 
representation will apply for [ I ] . 

Differentiating (B. 1-11) with respect to time 


[I] = 


(B. 1-14) 


but 


^(tr[K]) = /[x}^{x3dm) = / [x] ^ {x}dm 


+ f txl ’ [i 


[x3dm 


and recalling (B. 1-2) one may write 


[tx3^ 


;x] [x] dm = / 1 


{x]dm = I [x] Mx}dm = / tr [ [x] [x] ’'] dm =* / tr [ {x3 [x} ^Jdm 

M (B.1-15) 


Thus 


— (tr[]I]) - 2 / (x}Ux}dm 


(B.1-16) 


Similarly from (B. 1-12) 


~(tr[I]) = 2|(tr[I]) = 4 j {x]Ux}dm 


(B.1-17) 


Therefore (B. 1-14) may be rewritten in one of the following forms 


(B, 1-19) 





[I] = I It {tr[I])[l]-[I3 

or finally, 

[i] = i (A + B^- C)[l] - [i] (E.1-20) 

The values of [ I ] are readily available from any one of the above 
equations. 

B. 2 EKect of Coordinate System Translations on the Tensors of Inertia 

As implied by their definition, the values of [ I ] and [ I ] refer 
to a particular Cartesian coordinate system with origin at some specified 
point. 

The effects of a translation of the coordinate system on the tensors 
[ K ] and [ I ] and their traces will be investigated. In the first place, the 
formulation of the problem under the assumption that one of the coordinate 
systems is central is considered. 

B.2. 1 General Theorem of Parallel Planes 

This theorem relates the moments of inertia of a body about a plane 
through the center of mass to the moments of inertia about any other parallel 
plane, and viceversa- The proof given in this section is general and based 
on tensors of inertia. 

Suppose that the following values are known, 

M ^ total mass of the body in question 

[lo] 5 value of the tensor of inertia referred to a system (xq) 
centered at the CM ("central tensor of inertia") 
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} 


{6xo} = coordinates of the origin of the (final) system (x) after 

a parallel displacement with respect to the (Xo) frame 

From figure B.l one has immediately 

{x} =* {xo}-{6xo] (B.2-1) 

Thus, after substituting (B.2-1) in the definition of [I] given by equation 
(2.2-4), 


■I 


[I] = / [x}{xl^dm = I [xo}{xo}^dm - [6xq} 


y txo} ^ 


dm 


7' 


[6xo)^ /{xo}dm t6xo][6xo}^ / dm 
It is well known that 


/■ 


(B. 2-2) 


dm = M ; 


Therefore, finally 


Axol^dm = [O]^ ; jixo 

JtA 


}dm = [O] (B.2-3) 


[I] = tIo3 + M{6xo}{6xo}^ (B.2-4) 

The second- rank tensor {6xq3 {6xq 3 ^ will be called the "translation tensor" 
and represented by the (translation) matrix, 

[J^o] = [6xo3[6xo3' (B.2-5) 


The following theorem holds: The tensor of inertia with respect to 
any orthogonal frame parallel to the central coordinate system, is equal 
to the central tensor of inertia plus the product of the mass of the body and 
the translation tensor between the two systems. 

It should be noted that equation (B.2-4) is a generalization of the 
often-quoted theorem of parallel planes, according to which: The moment 
of inertia of a body with respect to any plane is equal to the moment of iner- 
tia with respect to the plane through the CM plus the product of the mass of 
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the body and the square of the distance between the two planes. 

IpL ~ H-opt Md p(^ 

From (B.2-4) it is obvious that 

(lol - [I] - M[^ol (B.2-6) 

which gives the central tensor of inertia as a function of the parallel tensor 
of inertia about any point and its translation matrix. 

Taking traces on both sides of the matrix equation (B. 2-4) 

trLii] = tr[Iol + Mtr[/iol (B.2-7) 

but (see fig. B.2-1) 

tr[/^o] = 6x^01 + 6x?oa + 6x®oa = d® (B.2-8) 

where 

d = distance between the origins of the (xo) and (x) systems. 

Thus (B.2-7) expresses what may be called the Theorem of Origin; 
The moment of inertia with respect to the origin of any system is equal to 
the moment of inertia with respect to the CM plus the product of the mass 
of the body and the square of the distance between the origin of the two sys- 
tems. The above can be applied to any point in E® and may also be obtained 
independently after substitution of (B.2-1) in (2.2-5). 

B.2.2 General Theorem of Parallel Axes 

As in the preceding section, the well-known theorem of parallel 
axes will be proved through a more general theorem relating the associate 
tensors of inertia. 

Substituting (B.1-12) into (B.2-7) 


M(5 


|tr[I] = |tr[Iol + Mtr[A''o] 
or 

tr[I] = tr[Io] + 2Mtr[A\o] (B.2-9) 

Recalling (B.1-13) and {B.2-9) one may write 

[I] = i(tr[Iol)[l] + M(trr 2^0 1)111 - [I] (B.2-10) 


From (B. 1-13) it is known that 

[TLo] = I (tr[Iol)tll - [lo] (B.2-11) 

Replacing the values of (B.2-10) and (B.2-11) in. (B.2-4) and after 
simplification 

[II = [lol + M[Ao] (B.2-12) 

v/here the "associate translation tensor" [ Ao ] is given by the matrix 

[Aol = (C6xb}M6xo3)[l] - l6xolt5xbr = [^][^]^ = -[6xp]'‘ 

(B.2-13) 


Hence the following theorem may be stated; The associate tensor of inertia 
with respect to any frame (x) parallel to a central coordinate system 
(Xq) is equal to the associate central tensor of inertia plus the product of 
the mass of the body and the associate translation tensor between the origins. 

Observe that as an application of equation (B.2-12) the moments of 
a body with respect to a given line or axis can be obtained. 

From (B.2-13) and fig. B.2.1, one may write, 




[ Aol = 


s 


— 6 X 01 6x 02 


d 


s 

2 


- 6Xoi 6 Xo3 

- 6x026x03 

d^3 


(B.2-14) 
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Thus from (B.2-12) the parallel axis rule follows immediately 

Iax = loAX + (B.2-15) 

The moment of inertia of a rigid body about a given axis through a 
point is equal to the moment of inertia about a parallel axis through the 
CM plus the product of the mass and the square of the distance between the 
two axes. 

Finally one should notice that equations (B.2-4) and {B.2-12) give 
exactly the same values for the products of inertia. This was expected by 
simple consideration of the definition of the quantities D, E and P. 

B.2. 3 Relation between TVo Non-Central "Parallel" Tensors of 
Inertia 

In this section the equation relating two non-central tensors of iner- 
tia [I] and [I'l referred to parallel systems (x) and {x') respectively, 
will be given. 

From (B.2-4) and the shift sign convention implied by fig. B.l, 

[I] = [lol + M{6xo][6xo3^ 

and 

[IM = [lol + M{ 6 xo]{ 6 xor 

Thus 

[I"l = [I] + M[{ 6 xa] [ 6 xol^ - t 6 xo}{ 6xo]^l (B.2-16) 

or 

[K'l = [I] + M [^ol-[Ao] (B.2-17) 

If the shifts [6x] between the (x) and (x') systems are known 
[6x] = [6 xo}-{ 6 xo] {6x0} = {6 x}+ {6x0} 

which can be substituted in (B. 2-lG) and gives 
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[I'] = [E]4 M[(6x}[8x}U [6x3 {6xo3'-i- [6xo3C6x3'] (B.2-18) 

A consequence of (B. 2-18) is that to obtain the relation between two 
non-central "parallel” tensors of inertia [ I ] and [ H^], the position of at 
least one of them with respect to the central (xq) system is required, in 
addition to knowing the shifts between the (x) and (x') systems. 

Similarly to (B.2-17) one has for the associate tensor of inertia, 

[I'l = [I] I M [A'o] -[Aol (B.2-19) 

and the same conclusion mentioned above also applies here. 


B. 3 Effect of System Rotations on Tensors of Inertia 


In this case one is seeking the tensor transformation 


[I] -^> [I] (B.3-1) 

where R is the orthogonal matrix of the inapping between the coordinate 
systems (x) and (x), to which the two tensors are referred, that is 


[x3 


R 


[x3 or [x3 = R[x3 


(B.3-2) 


Applying (B.3-2) to the definition of [ I ] 


[I] = f [x.) [xV dxn. = ^R [x3 [x3^R^dm = R f [x] [xY dm , 


which gives 

[I] = R[]I]R' 


(B. 3-3) 


This proves an important theorem in tensor analysis for transforming 
second-rank tensors under rotation. Observe the simplicity by which it was 
obtained here, using only the notions of matrix calculus. 


B.3. 1 Invariants under Rotation 


By definition, 
tr[n] = I [%y[x] dm 


Substituting (B. 3-2) above gives 
^[ 1 ] 


=jrcx] ^R^Rfx}dm = {x}dm = tr[I] 


Thus 


tr [ I ] = tr [ 1 3 
and recalling (B. 1-17) 

tr[T] = trin 


(B.3-4) 


(B. 3-5) 


Consequently the trace of the inertia matrix is invariant under rotation. This 
invariance of the traces of the inertia tensors makes an important check 
when equation (B.3-3) is numerically computed. 

K also follows that the determinant of the inertia matrix is invariant 
under rotation. Taking determinants in (B.3-3) and after consideration of 
the orthogonaliiy condition of R (R orthogonal ==> det R = det R^ = i- 1), 

det j = det R det [ I ] det R'f = det [I] (B.3-6) 

The same properties and transformations apply to the associate tensor of 
inertia. 

Equation (B.3-3) holds in general for any transformation under rota- 
tion of second-rank Cartesian tensors. Thus, for the associate tensor of 
inertia one can write. 


m = RriiR" 
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(B.3-7) 






The above relation may be proved independently after substituting (B.3-18) 
in (B.3-3) and considering (B.3-5). 

Recalling that the matrix R is orthogonal, from {B.3-3) and (B.3-7) 
it immediately follows that 

[I] = R^[i]R (B.3-8) 

and 

[I] = R'[T]R (B.3-9) 

It is Imown that R is the matrix of direction cosines of the three 
Cartesian axes of the (x) system with respect to the (x) system. 

B . 3 . 2 Applications 

(B lYansformations of variance-covariance matrices 

At this point one should notice that an immediate application of {B.3-3) 
to geodetic problems is the transformation of point variance- covariance 
matrices (VCM) between different Cartesian coordinate systems. In this 
case, if 

[Y} = R{X} 

then 

[Lvl = R[SxlR^ {B.3-10) 

where 

[Lx] - V CM with respect to the X system 

[ Ly ] = VCM with respect to the Y system 

When the transformation of variance-covariance matrices between Cartesian 
and curvilinear geodetic coordinates is intended, then the ibllowing linear- 
ized relation applies [Soler, 1976], 

{dxi, dxa, dXa) — (dX, d<p, dh) 


ini 


or 


dX' 


( dxi' 

d(p , 

- = jN 

dxa 

dh 

) 


,dxa^ 


(B.3-11) 


•where J is the Jacobian matrix of the transformation between Cartesian 
and curvilinear geodetic coordinates. 

Thus if the transformation between local VCM of Cartesian and curvi- 
linear geodetic coordinates is desired, one has 


where 


^^X<ph^ 




sinX 

cos X 

0 


(N + h)cos cp 

(N h h)cos <p 

J-^ = 

sincpcos X 
M+ h 

since sin X 
M + h 

cosep 
M+ h 


coscpcos X 

cos (psin X 

sin(p 

Note that in this case the second-rank tensor 



(B. 3-12) 


Sian" and the transformation matrix is not orthogonal, but of the Jacobian 
type. As a consequence, the following conversion of units will hold: 


(all m^) 


* 3 * 3 ^ 


(rad f 


{ rad f 
( rad f 


rad X m 
rad xm 

m2 

m 


Observe that the analytical form of equation (B.3-12) cannot be obtained 
directly through the classical propagation of errors approach, because the 
curvilinear geodetic coordinates are not expressable as simple functions 
of the Cartesian coordinates. 


(2) Momental and gyration ellipsoids 

A particular problem often mentioned in dynamics is the computation 
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of the moment of inertia with respect to a line of Icnown direction cosines 
{ 0£ } through the origin of coordinates . In this case the matrix R in (B. 3-7) 
will be reduced to a single row matrix [ a 

Thus, the moment of inertia with respect to the line (axis) when the 
direction cosines are known is given by 

Iax = [oj3"[I][od} (B.3-13) 

Similarly, the moment of inertia with regard to a plane whose 
normal has direction cosines {a}, is 

Ipu = [a3'[lUa} (B.3-14) 

H: one now considers the following relations 

{x] = p{oi3 (B.3-1G) 

where p is the radius vector of the {x} point, that is, 

? = (xrCx] = 

then 

[ot] = i tx} and [o:)^ = ^ {x}^ 

After substituting the above equation in (B.3-13), 

P^Iax= [xr[I][x] 

and taJting the value of the radius vector equal to 



one will have, 

[x3^[I][x] = M =? k = constant (B. 3-1(3) 

which is the equation of a quairic associate to the tensor [ I ] and the ixDlnt 
Ox; e is an arbitrary length introduced to keep the dimensions of equation 
(B. 3-16) coi'rect. 'iTiis quadric (known as the inei’tia or momental ellipsoid) 
has the properly' that the moment of inertia about tuiy radius vector is inversely 
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proportional to the square of that radius vector* 

Thus Iax will be maximum when p is minimum, and viceversa. 
Furthermore, if an, b„, Cn are the semiaxes of the momental ellipsoid 
at any point P with respect to which Ap, Bp and Cp are the principal mo- 
ments of inertia and if Ap < Bp < Cp , then ap > bn > Cn - 

The axes of this ellipsoid are the principal axes of the body at the 
point. The relation between the principal moments of inertia and the semi- 
axes of the momental ellipsoid is immediate. Using equation (B.3-1G) one 
may write 


or 

2 

ApX®i H- Bpx\+ CpX^a = k ^ 

Ap 

and finally, 

a®n = k/Ap , b\ = k/Bp , 

which gives 

Ap < Bp < Cp a„ > b„ > Cp 


2 2 

!l !l 

Bp Cp 


1 


c^ = k/Cp 


(B.3-17) 


(B. 3-18} 


If two principal moments at any point P are equal, the ellipsoid becomes 
a rotational ellipsoid; if the three principal moments are equal, the associate 
quadric to the tensor [TJ becomes a sphere, and every diameter is a princi- 
pal axis. 

It should be noticed that the longest and the shortest axes of the cen- 
tral momental ellipsoid (which is unique in a solid body) coincide in direc- 
tion with the longest and the shortest axes respectively of the material body. 
Thus the central momental ellipsoid resembles the general shape of the 
body, is elongated when the body is elongated and flattened when the body 
is flattened. 
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The reciprocal surface of the tnomental ellipsoid is termed "ellip- 
soid of gyration". Its equation is given by 


ifi + 2^+ ^ = i 

Ap Bp Cp M 


(B.3-19) 


and thus 

= Ap/M , b% = Bp/M , c% = Cp/M (B.3-20) 

Therefore 


Ap<Bp<Cp ==> a^<bs<Cg. 

Hence the (central) ellipsoid of gyration has the reciprocal shape of 
the body and its (central) momental ellipsoid. 


B.4 Final 'IVansFormation Equations 

The final transformation equations for the inertia tensors [ n ] and 
[ I ] after consideration of parallel shifts with respect to the CM and rotations 
of the Cartesian system, may be written by the following respective equations: 

[I] = R |^[Io] I- M[ AS] ] (B.4-1) 

and 

[I] = r|^[Io 1 + M[ Ao j] R^ (B.4-2) 

As an example of the utilily and application of the general equations , 
two particularly interesting problems may be reviewed: 

(i) Moment of inertia of a body of mass M about any straight line 
having equation 

-•= (B.4-3) 

CKi 0:3 Q!3 

with respect to a central system. Tlie above implies that the line has dlrec- 
tion cosines [a] and contains the j)oint (a, b, c); thus as a particular case 
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of (B.4-2) one may write: 




Iax = {a}'[[Iol+ M[ Aol] [O'} 


(B.4-4) 


or 


Iax = [o!}'[IoUa}+ M(0i}'[ Aol [a] 
where in this case 


[Aol - 


c® 


- ab 


a ' + c 


-a c 
-b c 
1 K' 


(B.4-5) 


(ii) Similarly, as an application of (B.4-1) the moment of inertia 
of a body about a plane containing the point (a, b, c ) and having normal 
vector of direction cosines [a] is: 


and 


UpL = ictV [[I 0 I+ M[ A\j 1 j [Oil 


(B.4-6) 


a 


[A\l = 


ab ac 

b^ be 


(B.4-7) 


The following subsections will try to clarify when a tensor of inertia 
is principal (i.e. , diagonal), depending on the orientation and origin of the 
system to which this tensor is referred. In other words , the condition under 
which the matrix equation 

r[[Io1 + M[ Ao1]r' = (B.4-8) 

is fulfilled will be investigated. Although not described here, the same pro- 
cess may be followed for the tensor [1 1. In conclusion it will be found that 
a body has infinite sets of principal axes, one at every point. However, a 
body has only one set of central principal axes. This intuitive assertion 
may be proved easily from (B.4-2) as follows: 


if) (5 


[I] = [lol only if 


[Aol = [0] ==> {6xol = {0] 
and 

R = [1] 


i.e. , the systems (x) and (xq) are coincident. 

Clearly the fulfillment of the same conditions will apply to the par- 
ticular case 

[•ij = rioj 


implying 

( Xp ) — ( X 0 p ) 

B.4. 1 Principal Axes Which Are Parallel to a Central System 


This section investigates the conditions under which a non-central 
tensor of inertia referred to a frame (x) parallel to a system centered at 
the CM is principal. 

By the parallelism assumption, one has U = [1]. It follows from 
equation (B.4-2) that only three cases are possible; 

(i) The central tensor of inertia is principal 

This is equivalent to solving the following diagonalfzation problem, 
[-10^+ M[ Aol = rij (B.4-t)) 

The above equation is true only when the matrix [ Aq ] is diagonal. That 
is, one of the three conditions given below must be fulfilled. 

6 xox = 6 xo 2 = 0 (B.4-10a) 

6 xoi = 6 xo 3 = 0 (B.4-10b) 

6xos = 6 xo3 = 0 (B.4-10C) 

Consequently the following theorem may bo stated; livery Cartesian system 
(x) parallel to the central principal axes (xoi.) and with origin in any of 
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its points is also principal- 

The following inequalities clearly hold from {B.4-9) 

Ap^Aop > Bp>Bop ) Cp^Cop 

That is, the minimum principal moments of inertia are obtained when the 
reference system has the CM as origin. 

Consequently, recalling (B.3-18) the momental ellipsoid of maximum 
volume is the central momental ellipsoid. 


(ii) Two central products of inertia are zero 

This condition implies that only one of the central axes is principal. 
Assuming, for example, that Do Eo = 0 (i.e. the x^axis is principal), 
one has 

r 

Aq “ Fq 0 

I 

+ M[Ao] ^ W (13.4-11) 


Bo 


0 

Co 


which immediately gives the following three conditions to be satisfied: 

Fo+ M 6 x 016 x 02 = 0 (B.4-12a) 

6 X 016 X 03 = 0 (B.3-12b) 

6 xoa 6 xo 3 = 0 (B.3-12c) 


Therefore , 


Fo 7 ^ 0 ==> 6 x 01 7 ^ 0 and 6 x 02 ^ 0 


The above implies that to fulfill the conditions b or c in equations ^^^4-12), 
6 xoa = 0. This means that the origin of the new (x p) system must be on 
the plane XoiXoa , that is, the plane through the CM and normal to the only 
central principal axis. 

Equation (B.4-12a) gives 


-F 


o 


iHK 


IVI 6xoi6xoa 




or 


OXqi 0Xo 2 - " jjj" 


(B.4-13) 


which is the equation of a rectangular hyperbola referred to its asymptotas 

Fo| 


M 


as coordinate axes and parameter K = 

Hence all the coordinate systems parallel to (xo) and with origin 
in the hyperbola given by (B.4~13), which lies In the x 01 X 02 plane, are 
principal (see Fig. B.2 ). Similar solutions can be obtained with 

F>o = Fo = 0 or Eo = Fo = 0. 


(iii) None of the central products of inertia are zero. 

In this case the given central tensor of inertia is not principal. Con- 
sequently one should solve the general diagonalization problem. 


[lo] M[Ao] = ru 

This implies the simultaneous fulfillment of the conditions 
Fo ^ M6 xoi6xo3 = 0 
Eo+ M6 xoi6xo3 = 0 
Do+ M6X026X03 = 0 

After solving for Sxoj (i = 1,2,3) above, one has 


(B.4-14) 


(B.4-15a) 

(B.4-15b) 

(B.4-15C) 


6x 


01 


- I ^0^’^ 
~ M M Do 


fiv - 1^0^° 


Eo 


(B.4-lGa) 

(B.4-lGb) 


_ IeoDq 


(B. 4 -IG 0 ) 


Therefore, in any body there are two specific points where one has a princi- 
pal system (Xp) which is parallel to any given centrjU frame (xo)^ The 
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two points are in opposite quadrants and equidistant from the center of mass, 
and their coordinates are given by equations {B.4-16). 

B.4.2 Eigentheory and the General Case with Rotations 

The most general case will be the computation of the principal tensor 
of inertia at a given point of a body when the central tensor of inertia is Icnown 
and no assumption of parallelism is made. 

This may be represented in matrix notation by the equation 

r-lj = R[T]R^ (B.4-17) 

where 

[T] = [Io] + M[Aol (B.4-18) 

Equation (B.4-17) expresses a well-known theorem in matrix diago- 
nalizati on theory, implying that any symmetric matrix [I] may be diagonal- 
ized by premultiplying and postmultiply ing it by some peculiar orthogonal 
matrix a’ld its transpose. 

It may be proved that 

R matrix of orthonormal eigenvectors of [ I ] 

[“IJ s diagonal matrix of the eigenvalues of [ I ] 

Thus, by the use of eigentheory the computation of ['IJ from [ I ] 
may be achieved. The problem is then reduced to finding the vectors 
{x} 1 (i = 1, 2, 3) such that 

[TUxlt = Xjx}i (B.4-19) 

Observe that equation (B.4-19) implies that the vectors [x}i (i = 1,2, 3) 
after being multiplied by the symmetric tensor [ I ] do not change their 
direction but only their magnitude, i-tu , any vector in a given direction, 
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after being premultiplied by a symmetric tensor, is stretched into a vector 
in the same direction but of proportionate size. 

This is not the case for tensors in general, which transform vectors 
into others of different magnitude and orientation. As is known, for example, 
the vectors transformed by orthogonal tensors keep their magnitude but change 
their orientation. 

It is well knoftn that by definition 

Xi (i = 1,2,3) “ eigenvalues of [I] 

[x3t (i = 1,2,3) eigenvectors of [ I ] corresponding to the 
X, eigenvalues 

In general, equation (B.4-19) is equivalent to 

[T (X)] {x} = {o3 (B.4-20) 

where the abbreviated notation 

[T(X)] = [T] - X[ll (B.4-21) 

has been used. 

The system of homogeneous linear equations given by (B.4-20) has 
a non-trivial solution [x} ^ {o} if and only if 

det[7(X)] = 0 (B.4-22) 

which is known as the characteristic equation of f I ( X) ] 

Equation (B.4-20) implies that the eigenvectors corrus] 3 onding to 
the different eigenvalues Xj are orthogonal to the row vectors of the ma- 
trices [I (Xi)]. Tlius the eigenvector direction cosines are proportional 
to the rows of [I (Xi)], i.e., to the rows of det [I (Xi)], consequently 
they will also be proportional to the minors of the elements of its rows. 

The left-hand side of (B.4-22) is called the characteristic polynomial 
fmd reprcaontod by (p( X). 
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Observe that 


(p(X) = det[I (X)l = X"" - aiX^+ aaX- ag = 0 (B.4-23) 

where the coefficients and the roots of tp(X) = 0 are related by the equa- 
tions 

ai = trlT] = tr[I] = X 1 +X 2 + Xg (B.4-24a) 

as = I ((tr[T]f - tr[Tf) - XiX 2 + X 2 X 3 + X 3 X 1 (B.4-24b) 

as = det[T] = det[n = XiXgXa {B.4-24c) 

Notice that starting with equation (B.4-19) one is referring to the eigenvalues 
and eigenvectors of the matrix [ T ] and not of H.( 1 ] as should be done 
according to equation (B.4-17), which expresses the initial diagonaliiiation 
problem . 

It will be shown now that the matrices 11 [ I] and [ I ] have identi- 

cal characteristic equations and therefore their eigenvalues and eigenvectors 
are the same. 

By simple consideration of the orthogonality properties of it, one 
may write 

det |^R[Tir^- X[ l]] = det[R[T]U^- RXlllR^] 

= det[R[[T] - X[ 1 ]]R'] 

= detRdet[^[Il - X[l] j detR^ 

= det[[Tl-X[ll] 

Thus, as a consequence, all coefficients in the characteristic equation are 
invariant under orthogonal transformation* They are called the principal 
invariants. Observe that the invariance under rotation of the coefficients 
tr [ I ] and det[ I ] was already proved independently in section B. 3. 1 . 

It is known that the values of Xi ( i = 1 , 2 , 3) for a real symmetric 


matrix are all real but not necessarily all distinct. 

Assuming now that Xi and Xa are two distinct eigenvalues of [ I ] 
and {x} and [y] their two corresponding eigenvectors, and making use 
of (B.4-19), one may write, 

[TUx] = XiCx] (B.4-25a) 

m ty} = XaCy} (B.4-25b) 

Premultiplying (B.4-25a) by {y]’^ and (B.4-25b) by {x}^ and subti’acting, 
[y}^[Il lx] - [x}^II] [y] = XilylUx] - Xa[x]Hy} 
but [ I ] being symmetric, the following equality holds 

ty}'m{x} = i;x}"[TUy} (B.4-2G) 

Therefore, finally 

(Xi - Xa) {x}^ [y] = 0 

By assumption Xi ^ Xs =p Xi - Xa 0 

thus {xl^{y3 = 0 "=c> [x] and (y] aro orthogonal. Hence the 

eigenvectors of a symmetric tensor corresponding to difl'eront eigenvalues 
are orthogonal (In general this is not valid for non-sym metric tensors). 

If the three eigenvalues are distinct, the three corresponding eigen- 
vectors form a unique orthogonal basis in if. The associated quadric to the 
[ I ] tensor will be an ellipsoid of tliree parameters. If two eigenvalues are 
equal, the associated quadric will be a rotational ellipsoid, implying that 
all the systems obtained by rotation of the axis with the only distinct eigen- 
value are possible. IE all the roots of the cubic equation (B.4-22) are equal, 
any system through the origin associated to the tensor are principal. 

As a practical application of the above eigentheory to the analysis 
of inertia tensors, suppose the general ease of finding the thi'eo principal 
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axes and moments of inertia of a body of mass M at a point {Sxqi, 6 xqs, 6x03) 
when the central tensor of inertia [Id is known. 

This is equivalent to determining the eigenvalues and eigenvectors 
of the matrix given by (B.4-18), namely 

[T] = [lol ^ M[ Ao] 


Therefore, the principal moments of inertia at (6x01 » 5xo2> 6x03) are the 
roots of the cubic equation 


Ao+M(6xo^ 2 '‘Oxo^al—X —Fo~M6xoi6xo2 

Bq+ M (6Xo~ 1 ** 6x0 3) “X 


s 


~Eo“M6xoi6xo3 
“Do” hi 6 xq 36X33 
Co+M(6xo^i+ 6 xo“3)“X 


0 


(B.4-27) 


Hence 



I 

H 

0 

0 

J 


Ap 0 0 


Xg 0 

= 

Bp 0 


s X3 


S 


(B.4-28) 


and the direction of the principal axes are the orthogonal eigenvectors corres- 
ponding to the values of Xj (i = 1,2,3). Notice that as mentioned above, not 
all the X 1 need necessarily be distinct. 


B.4.3 Some Geometric Interpretations 


(1) Foci of inertia 

K was proved in section B.4. 1 that all the Cartesian principal systems 
parallel to a given central system (with only one principal axis) liave their 
origin in a rectangular hyperbola, situated in the plane through the CM and 
normal to the only central principal axis. 


Clearly the system (xoi , Xq 2 » Xoa = Xopa) can be made coincident 
with the central principal frame (Xop) after a rotation 0 around the princi- 
pal axis Xop 3 (see Fig. B. 3). 

The rotated coordinate axes Xqi or Xo 2 will intersect the rectangu- 
lar hyperbola in two points symmetric with respect to the CM and with coordi- 
nates (6xopi, 0, 0) or (0, 6xops, 0 ) respectively. 

One may write the equation of the hyperbola in the new central princi- 
pal coordinate system (Xop) after considering the transformation of coordi- 
nates: 


txo) [Xo„] 

Thus 

{Xo^ = Ba'(0){xop] 

which gives 

Xoi = cos Gxopi - sinSxops 
Xo 2 = sin0Xopi + cos 0 xop 2 

Xq 3 ” X(jp3 


(B.4-29) 


(B.4-30a) 
(B. 4-30b) 
(B.4-30C) 


and similar equations hold for 6 xqi (i = 1,2,3). 6xos = 6xops = 0. 

Substituting the above values in the equation of the hyperbola (B. 4-13) 
and after considering tliat 


and 


■/ 


Fo = / XoiXosdm 


tan 20 


2sin6cos 0 
COS® 0 - sin* 0 


one has, after algebraic simplification, 

2M 6xo pi 6xo p 2 '• M tan2 0{6>^ o pi “6x®o ie)+ 9 


I T(x^; 


opi-x‘^op2)cim = 0 

{11.4-31) 






The above represents the equation of the rectangular hyperbola in 
the central principal coordinate system (xop). Observe that the variables 
are the coordinates (6xopi> 6xopa). 

By simple inspection of equation (B.4-31), if 




{x%pi-x^op2)dm >0 


I 


Opl 1 


>I 


Op23 


Bqp Aop 


(B.4-32) 


then 


5xopi = 0 


■which implies that the central principal axis will always intersect the hyper- 
bola along the axis of greater moment of inertia. 

Assuming (see Fig. B.3) 


Aop > Bop fixopa “ 0 (B.4-J3) 

then equation (B.4-31) gives 



(B.4-34) 

Clearly all the lines through the points ( 6xopi ,0,0) and (“6xopi ,0,0) 
in the plane Xopi xopg are principal at this point and with moments of iner- 
tia equal to Aqp* The momental ellipsoid at Ox is a rotational ellipsoid. 
Every body has six of these points known as foci of inertia, two in each prin- 
cipal plane. 


(2) Binet's theorem 

Assume the following diagonalization problem 

+m[Ao]Jh' = ru (b.4-35) 

where the diagonal elements of the central principal tensor of Inertia are 
given by 
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(B. 4 - 36 ) 


Aqp ^op WtS- Og * ^op Me Qg 

For simplicity the condition that the two equatorial moments arc the 
same is enforced ( aog = bog). 

The problem is then equivalent to determining the principal axes and 
moments of an inertia tensor [ I ] given by 


m = rioJ + M[Ao] 


(B.4-37) 


which is referred to a frame parallel to the central principal axes and with 
origin at the point P (6 xqi , 6xo2> 6x03). 

According to section B. 4 . 2 , this may be reduced to obtaining the 
eigenvalues and eigenvectors of the symmetric matrix [ 1 1 . Consequently 
the following equation in X must be solved; 


M(ao^g+ 6xo^a+ 6 xo^ 3 )“X -M 6xoi6xoa 

M (ao*g *■ Sxo'a'*" 6x0^ i)”X 

s 

Making the substitution, 

M { 6X(fi I- 6x0^ a >- 6x0 3) - X = MX 


-M 6 x 016 x 03 

-M 6X02 6x03 = 0 

M(Co*^ jj I 6x0 1 > 6xo‘ a)~X 

(B.4-38) 


(B.4-39) 


equation (B. 4 - 38 ) simplifies to 

n *■) 

-6xo"i t- (ao-g I X) -6x01 6x02 

-SxcTa'^ 


and the corresponding cubic equation becomes 


^6xqi 6x^3 
” 6 Xq3 6 Xq 3 

- 6 X 0^3 '■ (tio'g' X) 


= 0 


(B. 4 - 40 ) 


(ao“g+ly.)(Co^g+ X)- 6 xo\(Co^g+ X)- 6 xc?a(CoV ^)"Sxo%(ao"g+ X) = 0 (B.4-41) 
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Dividing through the above equation by (aog finally the following 

quadx'atic equation in X 

+ 5^0 3 _ ^ (B.4-42) 

ao^s+^ co^,*-"X. 

will have two roots "X l and 'X a- An ellipsoid of revolution will correspond 
to one of the "X, and a hyperboloid of revolution to the other. 

These surfaces will be confocal with 


6 xo\h 6 xq\ ^ fixp-’a ^ j 
adg 6 


(B.4-43) 


which is the rotational ellipsoid of gyration centered at the CM. As is known, 
the two confocal surfaces will intersect orthogonally . 

The orientation of the principal axes at the given point ( 6 xoi » 6 Xog, 6 xoa 
will be the same as that of the eigenvectors corresponding to the computed 
eigenvalues. It was shown previously in relation to equation (B.4-20) that 
for any Xi the direction cosines of the eigenvectors are proportional to the 
minors of the rows of the determinant given by (B.4-38), that is to say, to 

the determinant in (B.4-40). 

They are thus proportional to the quantities 

{aoV^)(Co\+X)-6xo%(ao^+'X)-6xo®a(coV^) ^ 6xoi 6xoa(cd’,+'X) ; 
6 xo^ 6 xo 3 {ao\^^lC}; dividing by (ap'^gi X){Co®g^ X) the above thr..e expressions 


reduce to 

6X03' 6 Xq2 

1 - - — 

cogt X aog' X 


dxp 1 6 xq 2 dxp X dxp 3 

Uq g ^ X Op g ^ X 


(Ii.4-44) 


but recalling (B.4-42), the above simplifies to 
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or finally, the eigenvector direction cosines are proportional to the quantities 


6xqi 6xq2 6xq3 

a© X ao X cq g+ X 


(B.4-45) 


which are the direction cosines of the normal to the, confocal surfaces passing 
through the point ( 6 xoi , 6 xoa » 6 xq 3 ) - 

Two of the principal axes at P are the normals to the respective mem- 
bers of the confocal surfaces which pass through the point. The third princi- 
pal axis must therefore be the axis mutually orthogonal to them, thus it will 
be on a plane through P parallel to the XoiXqs plane. I'wo of the principal 
moments of inertia are given by the eigenvalues Xi and Xa that may be 
determined from equation (B.4-39) after substitution of X. The third mo- 
ment can be computed using the parallel axis theorem as follows, 


M( 6 x^ 1 + 6 xo^ 2 "^ 3"^ 


This remarkable geometric finding that the principal inertia axes 
at any point P are the normals to the thi'ec orthogonal families of surfaces 
passing through P confocal with the ellipsoid of giration, was originally 
given by Binet [1813]. 

In the next section some important consequences of this seldom ref- 
erenced theorem in the diagonalization of symmetric tensors of second- 
order will be examined. 


B.5 Diagonalization by Rotation Matrices 

After the geometric conclusions derived from Binet's theorem, which 
unusually, links (dynamics and geometry, a novel approach (to the author’s 
best knowledge) for determining the pi-ineipal axes and momenta of inertia 
at any point of a body is presented hero. 




In [Soler, 1976] the theory for obtaining the rotation matrix R between 
the central and local frames {the axes of which are normal to the confocai 
surfaces) was studied. The problem here is thus reduced to computing 
the orthogonal matrix of the transformation between central and local frames 
from the available data. 

The initial values of M and ['IqJ will provide the ellipsoid of gyration 
parameters from which II may be determined at the siiecifiod point 

(6xoi 5 6 xq3 , Sxqs). 

Once R is known the matrix [ I ] given by (B.4-U7) can immediately 
be diagonalized 

R[i]R^ = rij 

By this simple geometric interpretation implicit in Binet's theorem, a dif- 
ferent method is found for obtaining the principal axes of inertia independently 
of the standard approach through eigenvalues and eigenvectors. 

A simple practical algorithm is given below for the particular case 
of (rotational) momental ellipsoid treated in the previous section. 


Central principal momonts 
Mass of inertia Point 


Given quantities: 


Central gyration 
ellipsoid para- 
meters 

Ellipsoidal para- 
meters 


Rotation matrix 


M 


App - Bpp Cqp ^ 


equation (B.4-36) 



'6Xqi , 6Xq3 , 6Xq3 ^ 

Transformation be- 
tween Cartesian and 
ellipsoidal coordi- 
nates [Heisltanen and 
Moritz, p. 228] 

l„ 

— ► (X, u ) 
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The value of R, in the case of rotational ellipsoidal coordinates was 
given in [Soler, 1976] as 




- sinX 


w„ 


sin/3 cosX 


u 




cos/3 cosX 


where 


cosX 


■■ — sin/3 sinX 

w. 


w„(u“±in^’' 


■cos/9sinX 


u cos/9 


w„(u®±IO^“ 


sin^ 


w. 


(B.5-1) 


/ u^±ifsin^/3 


, 1/2 


u ± 


(B.5-2) 


In the above equations the t- and - signs refer to oblate and prolate (rotational) 
ellipsoids respectively. In the latter case, E ='^co^g - acJg . The rows of 
are the direction cosines of the three principal (local) axes at the point 
P(6xoi, 6xoa, 6x03) " P(X, /8, iT). 

The principal moments at P arc the elements of the diagonal matrix 
obtained as a function of the given moments of inertia through the transforma- 
tion 

= m (B.5-3) 

where 

[I] = rio4-M[feof (B.5-4) 


The method can even be applied to the case of diagonalization of 3X3 symmetric 
positive definite non-singular matrices. 

Every symmetric matrix [ M ] may be expressed as the addition of 




i 


two other matrices, one diagonal [iVI^ and the other of the form [x} 

Thus 

[M] = + {x} [x]^ (B.5-5) 

To ['MIJ one may associate a central quadric that in the most general case 
will have three orthogonal families of eonfocal surfaces (equivalent to the 
ellipsoid of gyration in dynamics). 

At the point {x] the local moving frame will have the direction of the 
eigenvectors of [ M ]. Thus through this ijoint a eonfocal quadric to the one 
represented by will pass. (In dynamics this point will have principal 
axes parallel to the system associated to the tensor [ M ] ). 

llie point {x} then is found in form similar to the lo explained in 
section B.4.1. It will satisfy the condition 


and 


X1X2 

- 


X1X3 

= nii3 


XsXa 

= maa 


Xi 

^ /nHaniia 

(B. 5-6a) 

V niaa 


X2 


(B.5-6b) 



X3 

_ /niismsa 

(B. 5-6c) 

y iTiia 



Once the point {x] is Imown, the rotational matrix li between centi'al 
and local systems may be found as formerly described, using the parameters 
of the quadric associate to f'M'j. Finally, [ M 1 can be diagonalized using 
the stamlai'(! approaeli 
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R[M]R^ -- [-MJ 


The following subsection shows two simple examples where the rota- 
tion matrix R is immediately obtainable, although the method is general and 
may be applied to any 3X3 symmetric positive definite non-singular matrix. 

B.5.1 Examples 


(1) Spherical case 
Diagonalize the matrix [ M ] 



1 

2 

1 


1 

1 

2 


From equations (B.5-6) one has 

Xi = l ; X3 = l ; X3 = l 

Therefore the given matrix may be decomposed according to (B.5-5) as 


[M] = [ULi + [x] [x]' 


r "1 


- — 


— — 

to 


10 0 


111 



•< 



2 1 

- 

1 0 

V 

1 1 

1 

-1 


s 1 


5 1 


In this case it is clear that a sphere of radius 1 is associated to the matrix 
['Mij, ITius it will be easy to obtain the spherical coordinates of the point 
[x}. From Fig. R.4 it immediately follows 


sinX = 

^Ii 


cosX 



sini/) 



cosi/) - 


Vi 

vn 
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Fig, B. 4 , Local Frame ,C) and Base (01,02,63) 


Then the rotation matrix of the transfoianation between central and local 
systems is 


R= Ri(|-J/))R3(X + |) 


-sinX oosX 0 


-sinji)cosX -sini/)sinX cos^ 


cos^cosX cosi/tsinX sin(b 


_1 J_ 

nF ^j2 

_1 

" v/g y[G \j 

yj 


Thus the eigenvectors have the direction of the local base ( Oi ,ea ,63). F 
ally the given matrix f M ] is diagonalized as follows 


R[ M j [C 


0 0 

1 1) 


PMJ 




By the classical formalism of eigenvalues and eigenvectors, the same 


result is obtained. According to (B. 4 - 22 ) 


2-X 1 1 


det[M(X)] = 


s 


2-X 1 

2-X 


0 ==> {2-Xf - 3 ( 2 -X) + 2 = 0 


Therefore the eigenvalues are Xi = Xa = 1 and X3 = 4 . 

It may be shown that the unit vectors in the direction of the eigen- 
vectors corresponding to the eigenvalues X^ (i = 1 , 2 , 3 ) are the rows of the 
B matrix. 

Finally, it should be observed that the properties of the principal 
invariants described in section B. 4. 2 hold: 


(i) tr [ M ] = tr pM.] = Xi ^ Xa + X3 = 6 

(ii) ^ ^(trlMD^-trlMf j ■= ^ ^ (tr - tr pMJ^ 

Xi X3 t X3X3 + X3 Xi — 9 

(iii) det[ M ] = detpMJ = Xi Xa Xa =4 


( 2 ) Ellipsoidal (rotational) case 


Diagonalize the matrix [ M ] = 


2 1 
2 


0 

0 

2 


Following a procedure similar to example ( 1 ), one has 

Xi = Xa = 1 ; X3 = 0 

Thus 

[Ml 


1 

0 

“1 


r 

1 

1 

— 

0 


1 


\ 


1 

0 



2 * 




0 
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In this case the quadric associated to the matrix f’M L] is 


Xi+ xl+ 2x3 = 1 


P J 2 

xi t- Xa ^ X3 _ 


= 1 ilattened ellipsoid 


of revolution with parameters a = b = 1 and c = 




Thus the linear excentricity E is given by 

E =- 


E = -y a“ - 


c - — 


X[2 


The curvilinear ellipsoidal (rotational) coordinates at P(l, 1,0) may be obLainetl 
using equations (G-.S) in [lloisluinen and Moritz, 10G7, p. 228] 


M •» O ^ 1 /- 

u = (Xi > X5 I- X3 - E") ( 2 + gVJ 
_ x^u-+ fJ 


1 + 


4i«rx^ 


(xi > x’j. r - ir)*' 


3 

2 


tanjS = 


u^ 


- 0 ^ = 0 


f- Xa 


tan X - — = 1 ==5> ^ “ 7 

Xi 4 

Thus through the ))oint P(1,1,0) passes a confocal ellipsoid with parameters 




and E " — 

V2 




and from (B.5-2) it immediately follows tliat - -j - 

Thus the rotation matrix which makes the central (principal) system 
parallel to the local (moving) frame at P{Xi, = l,Xa - l,Xa - 0) P(X= ir/4, 
J3 = 0, u ~^js/yj2) can be obtained from equation (B.5-1). 


IA/2 1A/2 


IL 


0 


l.\/2 


1 

0 
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and finally 


lUMliC = 


0 0 

2 0 

3 


Observe that when the confocal ellipsoids are flattened, the resulting ele- 
ments in the diagonal (eigenvalues) are obtained in ascending order 


^ ^ Xg Xa 


and correspond b the directions of the tj, ^ and ^ local axes respectively. 

It may be shown as in example (1) tliat the principal invariant proper- 
ties hold. 

By the formal approach of eigenvalues imd eigenvectors the same re- 
sults are derived. 

Eigenvalues; 

1 

/ Xi - 1 


2-X 1 0 


2-X 0 

s 2-X 

- 0 (2-X)l(2-\f- 1 ] - 0 


Xa — 




Eigenvectors: 


A .1 = 1 : t M ] {x] = {x} ==5» I I ==> unit vector ( ^ ^ 


Xg = 2 : [ M 1 {x} = 2 [x] 


)xg = -Xi j 

Xl = Xg = 0 ( 


any Xa 


) 


^2 ,J2 

unit vector (0,0, 


{ Xl = Xg ) 

X 3 - 3 : [ M ] [x] = 3 [x] ==> I ^ _ Q I ==> unit vector ( ~ , — , 


V2 


In the general case the quadric associated to the matrix [ M 'j will be 
an ellipsoid of three parameters (i.e. the tliree elements in the diagonal 
will be different). Although the same approach may be followed, the intro- 
cluetion of general orthogonal ellipsoidal coordinates is I'equired and 



consequently the computation of R is more involved. 

By the same procedure explained above the errors ellipsoids from 
point variance-covariance matrices maybe obtained . Observe that the re- 
sult has strictly a geometrical meaning of difficult statistical interpretation. 
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APPENDIX C 

{Referenced in Sections 2.3.2, D. 1 and D. 1. 1) 

C. 1 Velocity and Acceleration in the Inertial System 

In this appendix the following notations will apply; 

(Xi.XgjXa) inertial (space fixed) system 
(xi,Xa,X 3 ) body fixed rotating system 

■ angular velocity components along the (X) system 
{ 6X } - coordinates of the CM with respect to the (X ) system 

It is well known that the transformation of coordinates between the 
two systems is given by the matrix relation 

[x} = r{x}+ lex} (c.i-i) 

where R is the rotation matrix between the (x) and (X) systems. 

C. 1. 1 Velocity 

The derivative of the vectors and matrices with respect to time will 

be denoted as usual by a dot. Let 

* 

[x] 1- velocity components of any point relative to the (X ) system. 
From (C.1-1) 

[X] = R [x] I- R lx] + {6X} (C. 1-2) 

but 

[x] ^ rUx} - RM6X] (C.1-3) 


ISO 


thus 


5 


[X] = rRT [X-6X] + r{x}+ [6X} 
where the "angular velocity matrix" is given by 
[OJ]x - HR" 

Thus finally 

[X] = [w]x[X-6X] t r[x}+ [6X} 
C.1.2 Acceleration 


Differentiating (B. 1-G) with respect to time, one has 
{X} = [w]x[X-6X}+ [m]x{x-6X}+ R{x}+ R{x] t {6X] 
or, after simplification, 

{X} = [wJx{X-6X]+ [£ofx[X-6X] h 2Rtx} ^ R {x} t [6X] 

but from (C.1-5) 

R - [(u]xR 
Thus finally, 

{X] = [a?lx[X-6X} i [wrUx-5X} (■ 2[w]xR[xl I R{x}i t6X} 
where the different terms in (C.l-G) and (C.1-8) have the following 

[WIxCX 6X} I j velocity ) of a particle P of the body in 

[i]xlX-5X)i I acceleration I the (X ) system resulting from 

( velocity v 

> of the (x) system 
acceleration ) 

[^]x[X-6X] ^ the centrifugal acceleration in the (X) system o 

P stationary in the (x) system 


(C. 1-4) 

(C. 1-5) 
(C.l-G) 

(C.1-7) 


(C. 1-8) 
meaning: 


■ a particle 
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2[w]xRtx} 


■v 

Coriolis acceleration of the moving particle 


R{x] = tX'} 


RCx}^tX'} 


, veloeito 

. 

ac(-valeration 


of a particle P in. the (X) system con- 
^ sidered relative to the (x) system 
) treated as stationary 


{fix} 

1 ( velocity 

{fix} 

) ( acceleration 


I of the CM in the (X) system 


C. 1.3 Velocity and Acceleration under the Assumption [SX] - 0 


If OX s Ox [6X] = {0} and [X] = R{x] , then equations (C. 1-6) 
and {C. 1-8) simplify to 

{X] = [toJxCx] + R{x] = [wlxtx] -1- [X'} (C.1-9) 

CX] = [wlxCX] + [WlKx]-!- 2[w]xR{x]+ R{x] 
or 

CX} = [W]xCX]^- [cojl[xl + 2[wlx{x0+ {XO (C.1-10) 

C.2 Angular Velocity Matrices 

The angular velocity matrix [w]x respect to the (X) system was given 
in (C. 1-5) by 

[coh = RR" (C.2-1) 

It may be proved easily that [W]x is a skew-symmetric matrix. Differen- 
tiating 

RR^ = [1] (C.2-2) 

one has 

RR^ + RlV = [ 0 ] = ^ RR‘ - - RR^ = -(RR^)^ 
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Thus 


[<^U = -Wl (C.2-3) 

The angular velocity matrix with respect to the body-fixed or rotat- 
ing system may be obtained after transforming the tensor [w]x under rota^ 
tibn (see section B. 3). 

K 

tx3-2^(x} or tX} = R{;x} (C.2-4) 

then 

[wlx = R[(0]R^ =5* [W] = R'^[^]xR (C.2-5) 

and substituting equation (C.2-1) above 

[w] = R^R (C.2-G) 

From (C.2-6) one immediately has 

R = R[C0] (C.2-7) 


which are the Poisson equations, and with 


R = 


Q£a «3 

Pz Pa 

7i 7z 7;i 


are generally given in the form 

dffit cloj 


^ ^ = odaWi - ttiOJa ; ^ = axCOg - 


da.. 


a _ 


and witb similar equations for V i = 1,2, 3 

In. particular if {a} are the direction cosines of a line respect to a moving 
fcame with angular velocity [io] one has 

[a]'^ « [ay[m] {C.2-8) 
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C. 3 Euler*s Kinematic (Geometric) Equations 


the so-called Euler's kinematic equations give the components 
(OijtOgjCOs) of the angular velocity as a function of Euler's angles tp, 0 
and^. The term "geometric equations" is applied by some authors. 

Assume that the transformation between the inertial and body-fixed 
systems is given by (see Fig. 2.3) 



11 

(C. 3-1) 

with 




H=R3(-<P)Ri(-0)Rs{-0) 

(C. 3-2) 

thus 




= R3(i/))Rx(0)Rs(cp) 

(C. 3-3) 


The components of {oj} along the (x) frame may be obtained from 
equation (C. 2-6) 

[CO] = ]RT IR (C.3-4) 


where 

* ■ 


K = Eg (-C0) Rx (-0) Ra (-0 ) -i- Ra (-<P) Ri ("8 ) R3 (-0 ) + R3 (~«p) Hi (- 9) Rs H ) 

(C.3-5) 

Considering that id general, for any angle x one has 


Ri(-X)Ri(x) = Ri(X)Rt(-x) « x[a]i Vi = 1,2,3 (G.3-6) 

where 

Rx(-X) = Rl(X) (C.3-7) 

and [a ]i are the antisymmetric matrices corresponding to the unit vectors 
which have the direction cosines of the three Cartesian axes. That is: 


= U, 0, 0}=^ [a lx 


0 

0 

0 


0 

0 

1 


0 

-1 

0 


(C. 3-8a) 


1B4 





i = 2=? {aF = {0, 1, 0} =^> [a]s = 


0 

0 

-1 


0 

0 

0 


i = 3=^ {a]T = {0, 0, 1} [a]3 = 


0 -1 
1 
0 


1 
0 
0 

0 

0 0 
0 0 


(G.3-8b) 


(C. 3-8c) 


From (C.3-6) it obviously follows 

rJ(x)Ri(h) = Ri{x)RI(x) = x[all Vi = 1,3,3 
Therefore from (C. 3-6) one may write 

Ri(-X)= X[aliRi{-X) Vi = l,2,3 (C.3-9) 

which can be substituted in (C. 3-5) giving 

il=(p[alalR+ 0Ra{"<p)[a]iRi(-e)R3(-0) + iLl3{-(P)^L{-B)[8L]aRa{-^) 

(C.3-10) 

and the angular velocity matrix takes the form 


— 

[ W] = ]R = Ra (5|)) Ri ( 0) ^ R 3 (cp ) [ala 1^ {"(p ) 9 [ ] Hi (~9 ) i [ 


HaH) 


Ri(0) 


= R3(>f}) 


0 -tp 0 

cp 0-0 
0 0 0 


Ri{-0)+ 0[al3 


HaH) 


AMtI SYMMETRIC 


-Q 3 COS 0-0 (p sin0 
0 -0 
0 




and finally 


[Wl = 


-(pooad-^ £p ain6cosi/)-0sitt^ 

0 -(p sin0sini|)-0coB0 

135 


{C.3-11) 


AHTUVMMFtniC 


0 


or 


Wi = ¥>311X03111^+ 0COS0 

■ • 

CO3 = (p sin0cos!/( - 8sini|) 

• ■ 

OO3 = <p cos 0 -h tp 

and in matrix form; 


(C.3-I2a) 
(C. 3-12b) 
(C.3-12C) 



sinSsini/) 

QOS ij) 

0 

1 

[ (P 1 

tto} = 

sin0coB0 

-sin0 

0 


0 


COS0 

0 

1 


^ 1 

By the same approach one 

may compute 

the compone 


(C. 3-13) 


inertial frame {X) 
[W]x = 


After matrix operation and simplification, one arrives at 



0 

coscp 

“1 

sin0 sincp 

1 

* \ 

Cw}x == 

0 

sin.(p 

-sin0 coscp 


0 


1 

0 

cos 0 


i 1 


(C.3-14) 


Although it is possible to obtain equations (C.3-13) and (C.3-14) by 
geometric considerations only, the approach followed here is general, and 
may be applied to any case where the rotation matrix H between the inertial 
and bod 3 '^“fixed frames is loiown. 
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APPENDIX D 


(Referettced in Sections 2. 3, 2.3. 2 and 2.4) 

Although recently many investigators have studied the mathematical 
theory of the rotational motion of deformable bodies about their center of 
mass, it is common to resort to the vector or tensor notation to formulate 
the principal equations. 

In this appendix the general Lagrange- Liouville equations extended 
to deformable bodies are derived only by means of matrix algebra, which 
one feels may also be used to advantage in formulating the basic equations. 

In another respect, the approach presented here is different from 
previously published ones. The point of departure is the fundamental axiom 
of mechanics called "Principle of balance and moment of momentum" that 
can be stated as: The time rate of change of moment of momentum about 
a fixed point is equal to the resultant moment about the point. 

One arrives at the same final results as previously published papers 
using varied approaches; i. e. [Kopal, 1968] bases his study on the fundamental 
equations of hydrodynamics, while [Tokis, 1974] starts from Cauchy’s first 
law of continuum mechanics; the analysis by Grafarend [1973 and 1977] 
originates in the so-called virial equations of second order. 

D. 1 Basic Dynamical Equations of Motion 

This section emphasizes primarily the general principles by which 
the differential equations describing the motion of an> material system of 
partic'les can be set up. Obviously the specification of the positions of the 
pai’ticles defines the body and its configuration. 
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The equations of motion of a body about its center of mass with respect 
to a space-fixed (inertial) system denoted by 



(Xi, Xs, X 3 ) = (X) 



is 

{l3x = {h3x 


(D. 1 - 1 ) 

where 

CH]< - 1 [H], 



and 

[h], =1 f [XlCx)dm 

■ column vector of total angular 
momentum about OX 

(D. 1-2) 


Cl3x = f [X]{fb}dm 
J M 

5 column vector of total body 
torques (about OX) acting on 
the body 

(D. 1-3) 


[ f ^ ^ bo^ force per unit mass exerted on each par- 

ticle of the body (see section D. 3-1) 

If = [ 0 } ^he system is said to be "isolated”. In this ease one 
has a torque- free body and the angular momentum vector is constant (i. e. , 
{H}x = 0 ==> {h]x = constant), having a fixed direction in inertial space, 
hence forming an invariant axis. 

Let a second rotating system (xqi ,Xo3,xoa) = (xo) be introduced, 
centered at the CM and not necessarily fixed to the body, but related to the 
inertial frame at any instant through the orthogonal rotation matrix R, i. e, , 

{xq} — > [X} or [X3 = R[xo3 (D. 1-4) 

(In what follows the vectors referred to the (xq) Cartesian system will be 
unsubscribed) . Then, 

[I-I3x = R{h3 (D.1-6) 

and 
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[H]x = R [H} + R{H} 


(D.1-6) 


Therefore, equation (D. 1-1) may be written. 

CL]x = R[H}+ R[H] (D.1-7) 

but 

[Llx = R{L] 

Thus 

[L] = RTR{h}+ R^R[ii} (D.1-8) 

The matrix 

R^R = [0] (D.1-9) 

is the angular velocity matrix with which the system (xo) is rotating. The 
components of (0] are referred to the (xd) frame. 

Observe tliat applying the theorem of section B.3 for transforming 
second-rank tensors under rotation, one may v/rite 

[_0]x = R[_0]R'^=RR^ {D.1-10) 

which gives the angular velocity matrix of the rotating system (Xq) with 
respect to the fixed frame (see also Appendix C). 

Substituting’ (D. 1-9) in (D. 1-8) and recalling the orthogonality con- 
dition of R, it follows 

{L] = [0][H]+ [H] (D.1-11) 

which is the equation of motion of the body with respect to the moving (rotat- 
ing) system when OX ^ OXq CM. 

hi the particular case when the rotating frame (xq) is body-fixed 

{0] = [o)} angular velocity vector of the rotating body (D. 1-12) 

and 

{T-3 = + [H] (D.1-13) 

this being the basic equation of motion in a rotating reference system tied to 
the body. 
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D. 1. 1 Total Angular Momentum 


The total angular momentum with respect to the fixed (X) system 
was defined in (D. 1-2) by 


[h)x = f [X]Cx}dm 

•/M 


The absolute velocity {X} may be expressed as (see Appendix C) 

{X3 = [W]x(X-6X] + R[x] + [6X] (D.1-14) 

Let us obtain the value of along the rotating frame (xo). Using the 
law for transforming second-rank tensors under rotation discussed above, 
one may write 

[X] = R[Xo]BJ 

Thus after multiplying the above matrix by equation (D. i-14) and referring 
the resulting free vector components to the (xq) system, one has 

[X][X3 = R[xo][CO]{xo] + R[3^][W]{6 xo3 + R[XplCxor3 ~ H[xolt6x3 (D. 1-15) 

Recalling equations (D- 1-2) and (D. 1-5) the value of the total angular momen- 
tum referred to tb.e (xq) system Is 


[H] = R’"Ch}x - .^R"[XUXldm 
*4i 


(D. 1-16) 


where [20 [X] is given by (D. 1-15). 

Therefore, after (D. 1-15) is substituted in (D- 1-16) the total angular 
momentum [H] may be given in the form 


[H] = Ch}t« + [h] i- [IIIrot 


(D. 1-17) 


when 


LOO 


(a) [H}th is the contribution to [H] due to the translation of the 


(D. 1-19) 


origin Oxq and is expressed by 

[hItr = / [Xol[feo]^dm{co} - y [xo]{6xo3dm (D.1-18) 

*'h *^h 

Clearly if {6 xq} = {O} Oxq ^ CM and 

[H]tr = {0} (D.1-19) 

(b) [h] is the relative angular momentum given by 

[h] = / [Xol[xor3dm (D. 1-20) 

Jn 

where {xor3 ' relative velocity of the body particles with respect 
to the (xo) rotating frame 

This contribution to {h} is due to the proper motion of the 
particles themselves besides the motion of the body as a whole 
due to rotation. 

(c) Finally [hJrot angular momentum caused by the rotation 

and has the form 


(D. 1-20) 


= f [Xol[Xo]^dmCco} = [Io]{oo} 


(D- 1-21) 


where by definition 


[lo] = 


/ [xoltxol’^dm 


is the central tensor of inertia associated to the system (xq). 


Hence, the total angular momentum respect to the rotating (Xo) system 


when Oxo s CM is given by 

[H] = [Io][£0}+ [h} 


(D. 1-22) 



D. 2 Classical Lagrange-LiouvUle Equations 


The equations of motioa of a changing body may be established by 
substituting equation (D. 1-22) in. (D. 1-11). The final differential matrix 
equation 

[L} = [Io][w} + [IoUw} + [0][Io]{6o}+[j9]{h]+ [h} {D.2-1) 

is not restricted to the rigid body case, as clearly can be seen by the pres- 
ence of the time derivative of the inertia tensor and the relative angular 
momentum vector [h]. Thus, equation (D.2-1) may be applied to more gen- 
eral situations, such as quasi-rigid bodies having particles moving among 
themselves, gyrostatic motions, expansion or contraction of systems of 
particles, etc. 

All the vector components in equation (D. 2-1) and the associate ten- 
sor of inertia refer to a central reference system of coordinates (xq) not 
necessarily attached to the body and rotating with angular velocity [6}. 

The vector {cjo} represents the angular velocity of the body about the axes 
of reference. Obviously when the axes of reference are body-fixed, [0] = {oo}. 

Although equation (D.2-1) is usually credited to Liouville [1858] 
because it appears that he formulated it explicitly for the first time, Tisserant 
f Mecanique Celeste , Vol. II, p. 500] remarlcs that the so-called Liouville 
equations are implicit in the general equations of motion of a dynamical sys- 
tem given already by Lagrange [1815]. 

m effect, the Lagrange equations may be expressed by 
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L . 

^ dtVSlOa, 


or in the abbreviated matrix notation used in this work, 


(D. 2-2) 


where 


and 


- 1 f[ini 


v/here 


{xj dm " total kinetic energy of the system (D. 2-3) 

« 

[x] " absolute velocity oE the particles. 

Assuming the origin of the (x) system at the CM {i.e. (x) (Xq) ) 
lx} = {Xq} = rWltXo} + {xbr} - [xor] 

{xbr} ^ components of relative velocity w. r. t. the (xo) frame, and 


Cxb}^[Xo] * [w3‘^[lipl[Xp]'' [co] + 2([to}^[Xol[xoJ) i- [xor}'^ [xof] 

Thus integrating the above equation over iiie mass of the body and making 
use of (D. 1-19) and (D. 1-22), 

T = I [ lol [to] + [w]^ [h] ‘‘‘ I t^or}^ {xor 3 dm 

where the value 

^ / [xoi,]’^ [xor } dm 5 kinetic energy of the relative motion 


(D. 2-4) 


Hence, taking the partials of equation (D.2-4) with respect to {to} 

[Iq]{co}+ Ch} (D.2-5) 

and finally computing the derivative respect to time of the above equation 
gives 
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(D. 2-6) 



m. 

{|^} = [iolCio}+ [IolCw}+ {h] 

After replacement of (D.2-5) and (D.2-6| In {D.2-2) one arrives at equation 
(D.2-1), which from now on will be referred to as the Lagrange-Liouville 
equation. The recent popularity of the Lagrange-Liouville type of equation 
in geophysical literature should be credited to [Munlc and MacDonald, 1960, 
p. 10] who applied them to the solution of some geodynamic problems. 

As mentioned previously, the Lagrange-Liouville equations {D.2-1) 
refer to a Cartesian system with origin at the CM. A more general expres- 
sion with respect to a parallel non-central system may be obtained by making 
use of equation (B.2-12), namely 

[I] = [Io] + M[Aol 

Consideration of the above will introduce the following contribution to equa- 
tion (D.2-1) 

M[il[Aoltw}+ M[Ao][w} + M[Aol[ob] (D.2-7) 

where clearly the velocities of the origin of the new frame are involved 

* 

through the matrix [ AqI • 

A final remark should be made here about the Lagrange-Liouville 
equations as given by (D, 2-1). They are functions of the associated central 
tensor of inertia [ Iq] • Recalling the relationship between the tensors [ I ] 
and [1] and their time derivatives as expressed by equations (B.1-11) and 
(B. 1-14), a new form of the Lagrange-Liouville equations in function of the 
tensor [ Iq] may he obtained easily, 

[L] = ^33 1^333) fw} “ [ ^ol fd)} '1* ( I 0X1+ Io2a'*' ^033 )[to} 

- + [ eiCb] + £h] 
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(D.2-8) 


D. 2. 1 Selection of Reference Coordinate ^sterns 


In tile application of the Iiagrange-Liouville differential expressions, 
various approaches are possible in order to select the reference coordinate 
system. It is a common practice to choose the reference system in such 
a way as to facilitate the maximum simplification of the equations. Several 
options follow. 

(i) Absolute body-fixed reference system 

For convenience it is usual to select a reference system fixed to 
the body in some way. It is known that in this particular case [6} “ {w} ; 
thus from (D.2-1) 

[L] = [Io]C<I!i3 + + taj][Io]{w} H- [COUhl + [h] (D.2-9) 

Observe that in this instance the tensor of inertia associated with the refer- 
ence coordinate system will vary if the mass configuration of the body changes. 
Thus the value of [ Iq] should be a known function of time. Hence the com- 
ponents of the tensor [ Iq] may be considered to be composed of a constant 
part [ I ]e and a contribution [ AT ] (time dependent) which will account for 
the differential variations in the tensor of inertia due to the changes involved 
in the body . Therefore, 

[lol - [I]e + CAT] (D.2-10) 

Similar reasoning may be applied to the relative angular momentum 
vector that should be expressed as 

[h] = Ch]E -I- {Ah} 

A note of precaution. Although in a nou-rigid body one may find such 
types of absolute body-fixed axes, the situation is complicated greatly in 
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tlie case of the real earth, primarily because the crust moves and any bench 
marks fixed to it defining the coordinate axes will inevitably change with 
time. 

(ii) Central principal axes as reference system 

I£ the instantaneous positions of the central principal axes are taken 
as the axes of reference, then Dq = Eo ~ Fo = 0 a-nd, the tensor of inertia 
win be diagonal. This simplifies considerably the explicit form of equations 
(D.2-1). Nevertheless, observe that the axis to which the diagonal tensor 
is referred will move if transport of mass is assumed. This is a major 
disadvantage of this option. 

Thus only at the initial Instant 

{0} = [m] (D.2-11) 

but for any subsequent time [6] = {co+ 6(o] where the vector {6a)) represents 
the variation in the components of {o)] as a consequence of the motion of the 
principal system itself; that is, the angular velocities with which the central 
principal axes (axes of figure) are separating from the Initial reference 
system. 

Hence equation (D.2-1) will reduce to: 

{L} = ^lo){(I))+ rio]{co}+ [a^Sa)]rioJ{w}+ [0)+6co]th}+ (h) (D.2-12) 

All the values in the above equation are referred to the instantaneous central 
principal axes. 

In this instance, as in (i) the principal moments of inertia must be 
known as a function of time, thus 

rioj = Me + rAy (D.2-13) 

In the particular case that the central principal axes move only differentially 
in the body, the following simplification applies 
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and 


[ATJrsIOI 


{D. 2- 14a) 


[6U)}fa{0} (D.2-I4b) 

(iii) Mean axes as reference system 

A third alternative otten mentioned, advocates the selection of the 
reference system in such a way that at the instant of selection 

[h] = [0] th) = [0} 

simpHy ing equation (D. 2-1) accordingly to 

[L} = [lo3[oj3+ [IolCw}+ [eiEIolCt*)} (D.2-15) 

Although the imposition of the above condition will "freeze" rotations (currents) 
of the moving particles, the body still can have radial deformation (see 2.5.6). 
This set of axes has been called "mean axes" [Tisserand, 1891, Ch. 30] but they 
have the disadvantage that they are not uniquely determined [Routh, II, 1905, p. 21]. 


D. 3 Lagrange- LLouville Equations Extended to Deformable Bodies 

In section D. 1 the total torque acting on the body was defined by 
equation (D. 1-3) 

[l3 = 

where {fi>} represents the body forces per unit mass exerted on each particle 
of the system. 

These are not the only forces that may act on the deformable body. 

It is also possible to consider surface forces per unit area {f^} which will 
create surface torques of the iype 
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(D.3-1) 


M - 

Body and surface forces are the resultant effect of several components origi- 
nated in different ways (i. e. , directly applied forces, reaction forces due 
to constraints on the system, etc. ) and which finally may produce deforma- 
tions on the body. 

The individual description of these forces and the mathematical formu- 
lation to compute them follow. 

D.3.1 Body Forces 

The body forces acting upon a body may be divided into two large 

groups: 

(i) External forces resulting fro-n the interactions with objects out- 
side the system (i.e. , gravitational attractions, cdectrostatic 
forces, etc.) 

(ii) Internal forces corresponding to different interaction of particles 
located in the interior of the body (i. e. , Newtonian gravitational 
forces, centrifugal forces associated with rotation, etc.) 

Body forces are expressed in general by the vector equation 

[f.l " (D.3-2) 

is a vector-operator and the scalar represents 
the total potential acting on the body. 

The value of W is given by the sum of all the contributing potentials 
in the particular problem at hand. 

For example, 


•where the symbol \ 


I. 


[x] {f ,} ds 
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where 


*= Vg + Vm + Vq + + Vjj + other potentials 

V E = potential of the stationary (equilibrium) state of undistorted 
configuration 

Vh = the total potential arising from the mass 

, Vj , Vg s gravitational potentials due to the attraction of the 
sun, moon and planets. 

D. 3- 2 Surface Forces 


Surface (contact) forces, also referred to as surface stresses, de- 
rive from the action of one portion of the material on another through the 
bounding surface (i. e. , forces from hydrostatic pressure, wind drag, fric- 
tion, viscous force, etc.) 

These forces are given by 


[fj = [orlW 


(D. 3-3) 


where 


[a] - 


Oil CTl2 CTl3 
CTga Og3 

3 


“ stress tensor 


(D. 3-4) 


[n} =5 vector normal to the element of surface ds and with 
magnitude equal to the area 

The corresponding torque, after consideration of the Gauss integral 
equation^ may be written 


The stress tensor of an isotropic elastic medium (one whose elastic pro- 
perties at any point are independent of direction) is given by 


tCJl = X(tr[e])[ll + 2/X[D] 


{D. 3-6) 


where the following notations are introduced 



®11 


■ 

5xi 

3xi 

■ 

3xg 

3xi 

• 

3xa 

3xi 

[e] = 

GgS ®23 



■ 

Bxg 

Bxg 

• 

3x3 

3xg 


5 ©33 


3 


• 

3X3 

3xa 


strain 

tensor 


(D. 3-7) 


tr [ e ] = 6o = dilatation (increase in volume per unit of volunve of the 
medium 


[ D] = 1(11^ I [xf + + [ ef )- deformation rate 

tensor (D, 3-8) 


X,}X ^ Lame’s elastic constants. These quantities are the most con- 
venient parameters to use in theoretical work. In applications 
to the earth they are methodically substituted by Love’s num- 
bers. Nevertheless, Lame's constants have the advantage of 
being directly related to other common elastic parameters, 
such as Young's modulus and Poisson’s ratio. 

In the case of a linear viscous fluid, one has 

[O] = -(p- Xtr[e])[ll i- 2pi[D] (D.3-9) 

where now 

p = hydrostatic pressure 
/i s coefficient of viscosity 

it can bo proved tliat in this Instanco 
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(D. 3-11) 


and equation (D. 3-9) may l)e written in the form 

[cr] = -(P+ |Mtr[e])[l]+ 2(i[D] 

In the particular case of perfect fluids fi = 0, thus 

[cr]=-p[ll (D.3-12) 

When applying the above equations , it is typical to begin with the assumption 
that the body initially is in hydrostatic equilibrium 

[c7l = -po[l] 

and then examine small deviations from these equilibrium conditions. 


D. 3. 3 Final Eqiiations 


The final Lagrange- Liouville equations applicable to a deformable 
body, when elasticity, viscosity, etc. are taken into consideration, is given by 


{l 3 + iLj = [To] to)} + [ io] [io] + [ 0'K lol -i- [01 [h] + {li} (D. 3-13) 


Ihe above matrix differential equation is in a form sufficiently general to 
impose few restrictions on deformation when the contributions to [Iq] and 
th] are properly computed. 

The torques {L] and [Lg] must be obtained for the particular assump- 
tions involved, and in general are expressed by 


tL] = 


= / 


■¥dm 


(D. 3-14) 


where V represents the total ijotential created by the bociy forces acting on 
the system, and 




CL.} = 




dv 


{D.3-15) 


wMeh. gives the equations of the torques governing the appropriate stress- 
strain relations. 

All the vectors and second-order tensors in equation {D. 3-13) refer 
to a pre-selected reference system with the simplifications described in 
section D. 2. 1 possible. 
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